{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 回帰分析"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.106460Z",
     "start_time": "2018-08-24T07:49:05.600877Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import stats\n",
    "import statsmodels.formula.api as smf\n",
    "\n",
    "%precision 3\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.119850Z",
     "start_time": "2018-08-24T07:49:06.107608Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "20\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>小テスト</th>\n",
       "      <th>期末テスト</th>\n",
       "      <th>睡眠時間</th>\n",
       "      <th>通学方法</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4.2</td>\n",
       "      <td>67</td>\n",
       "      <td>7.2</td>\n",
       "      <td>バス</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>7.2</td>\n",
       "      <td>71</td>\n",
       "      <td>7.9</td>\n",
       "      <td>自転車</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.0</td>\n",
       "      <td>19</td>\n",
       "      <td>5.3</td>\n",
       "      <td>バス</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3.0</td>\n",
       "      <td>35</td>\n",
       "      <td>6.8</td>\n",
       "      <td>徒歩</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>1.5</td>\n",
       "      <td>35</td>\n",
       "      <td>7.5</td>\n",
       "      <td>徒歩</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   小テスト  期末テスト  睡眠時間 通学方法\n",
       "0   4.2     67   7.2   バス\n",
       "1   7.2     71   7.9  自転車\n",
       "2   0.0     19   5.3   バス\n",
       "3   3.0     35   6.8   徒歩\n",
       "4   1.5     35   7.5   徒歩"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.read_csv('../data/ch12_scores_reg.csv')\n",
    "n = len(df)\n",
    "print(n)\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 単回帰モデル"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.124160Z",
     "start_time": "2018-08-24T07:49:06.121009Z"
    }
   },
   "outputs": [],
   "source": [
    "x = np.array(df['小テスト'])\n",
    "y = np.array(df['期末テスト'])\n",
    "p = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.233613Z",
     "start_time": "2018-08-24T07:49:06.125769Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAl4AAAF1CAYAAAA5ouTuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xt0VNlh5/vv1lsgkHiDBI2gAfEUCAn0KtH0C+gHjbrp7iSexHbidk+cm1znxmGWe5LppLMmsxwTJ3Ov1207ndhjJ44zjseY6/hemzjjbs+p0hMhhEAgnhIg8RRIvErP2vcPCdxN8xLUqaOq+n3W6oV0uqrOD20h/dapffY21lpERERExH0JXgcQERERiRcqXiIiIiIRouIlIiIiEiEqXiIiIiIRouIlIiIiEiEqXiIiIiIRouIlIiIiEiEqXiIiIiIRouIlIiIiEiEqXiIiIiIRkuR1gDuZOnWqzc3Njeg5+/v7SUlJieg5Jfw0jrFB4xgbNI6xQeN4fw0NDRettdMe5LFjsnjl5uaye/fuiJ6zra2NSJc9CT+NY2zQOMYGjWNs0DjenzGm/UEfq7caRURERCJExUtEREQkQlS8RERERCLElTlexph04O+BacAB4HeB94BC4MfW2rdH+5oDAwOcPn2a3t7esGa9aXBwkIMHD7ry2tEkLS2N2bNnk5yc7HUUERGRmOPW5PrfAH5qrf2GMeYrQBmQBhQBncaYr1trO0fzgqdPn2bChAnk5uZijAl74L6+PlJTU8P+utHEWktXVxenT59m3rx5XscRERGJOW691ZgAzBj5eApQALwP5AP+kc9Hpbe3lylTprhSumSYMYYpU6a4dlVRREQk3rlVvL4NLDTG7AKuApnARWA68MHI56Om0uU+fY1FRETc49Zbjf3AG9baIWPMmwy/zZgKtDA8z+vE7U8YedybADk5ObS1tX3k/w8ODtLX1+dS3OHXv5++vj4++9nPcurUKSZPnsw3v/lNPvnJT3LlyhVSU1P51re+xcyZMz/2vF27dvGlL30JgJ6eHsrLy/nqV7/KO++8w09/+lPy8/P5+te/fqv0vPHGG/zd3/3dref/wR/8AXV1deTn5/Puu+/eMdtPfvITtm3bxrRpw+u3/cu//AuO43zsWEZGxgN9LW7/+keLrq4uryNIGGgcY4PGMTZoHMPLreK1HngW+CLwIrAd2GKt/b4xZh3DV8Q+wlr7HsMT8CkqKvrYyvUHDx50fQ7W/V7/+9//PoWFhXzve9/jC1/4Aj/96U/58pe/TH5+Pn/+539OTU0Nv/Irv/Kx57300ku89NJLALz99tuUlpbS1dVFIBBgz549FBUVsXfvXkpKSgBITEy8lSUQCDB58mTq6+v57Gc/y9GjR1m2bNnHztHX18cf/dEf8elPf/qexx5EUlJSVC+WF83Z5Zc0jrFB4xgbNI7h49Zbje8Dk40x1cBhhud1ZRljaoH60U6sHytWrlzJJz7xCWC4nCxYsIDU1FQKCwv52c9+xiuvvHLP51tr+dnPfsaGDRvYs2cPTzzxBOfOnSMvL4+GhoY7PudnP/sZAwMDPPXUU2RnZ9+xdAFcu3aNb37zm5SVlfHGG2/c9VhbWxtPPvkk58+fp7S09IGu9ImIiEh4uHLFy1obYuRtww95I1yv/9Of/pSzZ8+G6+UAmDZtGi+88MI9H7NixQoAdu7cSU9PD2VlZQA0NDTwzjvv8MMf/pDXX3/9rs//xS9+QVlZGYmJifT09DB16lSampr41V/9VZqbm/na177G9773PQ4dOsT69ev5/Oc/z9mzZ0lLS+PnP/85r7/+OvX19axZs+Zjrz137lw+//nPs3XrVj71qU/hOM4dj1VUVPDkk0+yefNm/vIv/5KkpDG5a5SIiEhM0gKqo/Sd73wHv9/P1772NU6cOMGFCxeA4bcT33///Xs+95/+6Z/YunUrAJmZmSQkJDBz5kyuXbtGZmYmn/vc5/jggw/YtGkTH3zwAS+//DIZGRmsX78egCeeeIKWlpY7vvYzzzxz67WLi4s5duzYHY8BPPnkkxw5coRVq1Y98tdDREREHlxUXu7YtGlT2F/zQSbut7e384Mf/IAf/vCHANTW1nLw4EHeeecd9uzZc8+1r0KhEI7j3Jocv3r1ar7zne/wO7/zO7z77rt3nYdVUlJCbW0tlZWV7Nmzh89+9rN3fNxf/MVfsHz5cl544QVqamr43Oc+d8dj1lr++I//mHfeeYc//dM/5Stf+cp9/94iIiLRYGdjB9t3tdLZHSQ7K51tG/OoLMjxOtZH6IrXKHzjG9+gubkZn8+Hz+ejr6+PEydOsG7dOn784x/z27/923d9bnV1Nfn5+SQmJgKQnZ3NokWLKCkpYWhoiOLi4js+r7KyktOnT1NWVkZWVtatCfi3+63f+i3+6q/+iqKiIjIzMyktLb3jsW984xusXbuW3/u936O6upoDBw48+hdGRETEYzsbO3hrRzMd3UEs0NEd5K0dzexs7PA62kcYa63XGT6mqKjI7t69+yPHDh48yJIlS1w7ZzhXrk9KSrq11IQbfvGLX/D0009z7Ngx5s6dG/bXd/tr7aa2tjbdfRMDNI6xQeMYG6JlHMu/9HM6uoMfO56TlU7gi0+5em5jTIO1tuhBHhuVbzWOdXv37mX27Nmuvf7N5Seys7NdO4eIiEg06bxD6brXca+oeLlg+fLlrr7++PHjXT+HiIhINMnOSr/jFa/srHQP0tyd5niJiIhI1Nu2MY/05MSPHEtPTmTbxjyPEt1ZVF3xstZqL0GXjcU5fyIiIvdz8+7FsX5XY9QUr7S0NLq6upgyZYrKl0ustXR1dZGWluZ1FBERkVGrLMgZc0XrdlFTvGbPns3p06dvLVgaboODg1rFneGC6+aNASIiIvEsappGcnLyPRcofVTRcrusiIiIRC9NrhcRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERkQhR8RIRERGJEBUvERERiRk3btzg5MmTXse4qySvA4iIiIg8qitXrlBVVcWePXtIS0vj93//90lIGHvXl1S8REREJGpdunSJQCBAU1MToVCI/Px8ysvLx2TpApeKlzEmFfgmMBe4Cnwa+M9AIfBja+3bbpxXRERE4sP58+fx+/3s37+fhIQEVq1aRXl5OZMmTfI62j25dcXrGeCUtfbfGWP+FHgTSAOKgE5jzNettZ0unVtERERiVGdnJ47jcOjQIZKTkykpKaG0tJQJEyZ4He2BuFW82oAfGGNeBi4D3wHeB/IBP1AAqHiJiIjIA2lvb8dxHI4dO0ZaWhrr1q2juLiYcePGeR1tVNwqXmeAl621PzHGvAVMAU4C04EPgEyXzisiIiIxwlrL0aNH8fv9nDx5kvHjx/P000+zZs0aUlNTvY73UNwqXtuA7418vBNoAV4f+bMQOHH7E4wxbzL8liQ5OTm0tbW5FO3Ourq6Ino+cYfGMTZoHGODxjE2eDGO1lpOnjzJvn37uHTpEuPHj2ft2rUsXLiQpKQkzpw5E/FM4eLmXY1lwF6gFHgbKLbWft8Ysw749u0Ptta+B7wHUFRUZHNzc12MdmdenFPCT+MYGzSOsUHjGBsiNY5DQ0Ps378fv9/PxYsXmTx5Mi+99BL5+fkkJiZGJIPb3Cpefw181xjzG8BZhu9q/IoxphbYpYn1IiIictPg4CB79+4lEAjQ3d3NjBkz2Lp1K0uXLh2zy0I8LFeKl7X2PMN3Nn7YG26cS0RERKJTf38/u3fvprq6mmvXrpGTk8OmTZtYtGgRxhiv47lCC6iKiIhIRAWDQerq6qitrSUYDDJv3jxeeeUVcnNzY7Zw3aTiJSIiIhFx7do1ampqqK+vp7+/n0WLFlFRUcHs2bO9jhYxKl4iIiLiqp6eHgKBAI2NjQwNDbFs2TJ8Ph8zZszwOlrEqXiJiIiIK7q6uvD7/ezbtw+A/Px8fD4fU6ZM8TiZd1S8REREJKzOnTuH4zi0tLSQmJhIUVERZWVlZGZq/XQVLxEREQmL06dP4zgOhw8fJiUlhbKyMkpKSsjIyPA62pih4iUiInFjZ2MH23e10tkdJDsrnW0b86gsyPE6VlSz1nLixAn8fj8nTpwgPT2d9evXs3btWtLT072ON+aoeImISFzY2djBWzuaCQ4MAdDRHeStHc0AKl8PwVrL4cOHcRyHjo4OMjIy2LBhA4WFhaSkpHgdb8xS8RIRkbiwfVfrrdJ1U3BgiO27WlW8RiEUCtHS0oLf7+fcuXNkZWXxwgsvsGrVKpKSVCvuR18hERGJC53dwVEdl48aGhpi3759+P1+Ll26xNSpU6msrGT58uUxs49iJKh4iYhIXMjOSqfjDiUrO0vzkO5lcHCQ2tpaqqqquHLlCrNmzeK1115jyZIlMb/KvBtUvEREJC5s25j3kTleAOnJiWzbmOdhqrGrr6+P+vp6AoEAvb29PPbYY2zevJnHH39chesRqHiJiEhcuDmPS3c13tuNGzeora2lrq6O3t5esrOz2bBhA3PnzvU6WkxQ8RIRkbhRWZCjonUXV69epaqqioaGBgYGBli8eDEVFRX09/erdIWRipeIiEgcu3z5MoFAgL179xIKhVixYgXl5eVMnz4dgLa2Nm8DxhgVLxERkTh04cIF/H4/zc3NJCQksHLlSnw+H5MmTfI6WkxT8RIREYkjZ86cwXEcDh48SHJyMsXFxZSWljJx4kSvo8UFFS8REZE4cPLkSRzH4ejRo6SmplJRUUFxcTHjx4/3OlpcUfESERGJUdZajh07huM4nDx5knHjxvHUU0+xZs0a0tLSvI4Xl1S8REREYoy1lkOHDuH3++ns7GTixIls2rSJ1atXk5yc7HW8uKbiJSIiEiNCoRD79+/H7/dz4cIFJk2axObNm1m5cqW29RkjVLxERESi3ODgIHv37iUQCNDd3c306dN55ZVXWLZsGQkJCV7Hkw9R8RIREYlS/f39NDQ0UF1dzdWrV8nJyWHTpk0sWrRI2/qMUSpeIiIiUaa3t5e6ujpqamoIBoPk5uZSWVnJvHnzVLjGOBUvERGRKHH9+nWqq6upr6+nv7+fRYsW4fP5mDNnjtfR5AGpeImIiIxxPT09VFVVsWfPHgYHB1m2bBk+n4+ZM2d6HU1GScVLRERkjOrq6iIQCNDU1ARAfn4+Pp+PKVOmeJxMHpaKl4iIyBhz7tw5/H4/Bw4cICEhgcLCQsrKysjKyvI6mjwiFS8REZEx4vTp0ziOw+HDh0lJSaG0tJTS0lIyMjK8jiZhouIlIiLiIWstbW1tOI7DiRMnSEtL44knnqC4uJj09HSv40mYqXiJiIh4wFrLkSNHcByH06dPk5GRwTPPPENRURGpqalexxOXqHiJiIi4aGdjB9t3tdLZHSQ7K50/3LCQhSk9+P1+zp07R2ZmJs8//zwFBQUkJenXcqzTCIuIiLhkZ2MHb+1oJjgwhCFE+tVT+H9Uzz7Ty5QpU9iyZQsrVqzQPopxRMVLRETEJdt3tdI/MMDixIusSDpLRkI/XaF0mpIX8z9+5zXtoxiHXClexpjfAn5r5NOZwN8CC4BC4MfW2rfdOK+IiMhY0dfXx+Rrx/GlnWOcGeTc0Hiq+x7jdCgT02dUuuKUK8XLWvtN4JsAxpi/B/4XsBwoAjqNMV+31na6cW4REZFwuX1+1raNeVQW5NzzOTdu3KC2tpa6ujqKknvpGJrIB4OzOBfKAIb3UczO0t2K8crVtxqNMRnAHIavdL0P5AN+oABQ8RIRkTHrw/OzADq6g7y1oxngjuXr6tWrVFdXs3v3bgYGBli8eDFD0/L47++fJRgauvW49OREtm3Mi8xfQsYct+d4vQrsADKBk8B04IORz0VERMas7btab5Wum4IDQ2zf1fqR4tXd3U0gEKCxsZFQKMTy5cvx+XxMnz4dgHGTRn/VTGKX28XrV4HPAC8DqUALw1e/Ttz+QGPMm8CbADk5ObS1tbkc7aO6uroiej5xh8YxNmgcY0O0j2Nnd/Cux9va2ujp6aG5uZnjx49jjGHBggUsW7aMiRMncuPGjVu/x1ZNgn/81fkfeoWBiP+OexTRPo5jjWvFyxiTCUyx1nYYYxqArdba7xtj1gHfvv3x1tr3gPcAioqKbG5urlvR7sqLc0r4aRxjg8YxNkTzOGZnHafjDuUrb+IQdXV1HDx4kOTkZIqLiyktLWXixIkepIyMaB7HscbNK14vAv9z5OMa4DPGmFpglybWi4jIWLdtY95H5nhNT7hKQfJZsgd6OH48lYqKCoqLixk/frzHSSWauFa8rLX/+KGPLfCGW+cSEREJt8qCHKy1/Lef1pLT28bMxGskpqTxhO8p1qxZQ1pamtcRJQppAVUREZHbWGtpbW3lYoND4UAnE7ImUFa2kcLCQpKTk72OJ1FMxUtERGREKBRi//79+P1+Lly4wKRJk3jxxRdZuXKl9lGUsNB3kYiIxL3BwUGampoIBAJcvnyZadOm8corr7Bs2TKtMC9hpeIlIiJxq7+/nz179lBVVcXVq1fJzs5mw4YN5OXlYYzxOp7EIBUvERGJO729vdTV1VFbW8uNGzeYO3cuW7ZsYf78+Spc4ioVLxGJSg+zh57I9evXqampob6+nr6+PhYuXIjP5+Oxxx7zOprECRUvEYk6o91DT+TKlStUVVXR0NDA4OAgS5cuxefzMWvWLK+jSZxR8RKRqPOge+iJXLp0Cb/fT1NTE9Za8vPz8fl8TJ061etoEqdUvEQk6txrDz0RgPPnz+P3+9m/fz8JCQmsXr2a8vJysrKyvI4mcU7FS0SiTnZW+h330MvOSvcgjYwlHR0dOI5Da2srKSkplJaWUlJSwoQJE7yOJgKoeIlIFLp9Dz2A9OREtm3M8zCVeMVaS3t7O47jcPz4cdLS0njiiScoLi4mPV1lXMYWFS8RiTo353Hprsb4Zq3l6NGjOI7DqVOnGD9+PM888wxFRUWkpqZ6HU/kjlS8RCQqVRbkqGjFqVAoxMGDB/H7/Zw9e5bMzEyef/55Vq1apX0UZcxT8RIRkagwNDREc3Mzfr+frq4upkyZwpYtW1ixYgWJiYlexxN5ICpeIiIypg0MDLB3714CgQA9PT3MmDGDV199lSVLlmgfRYk6Kl4iIjIm9fX1sXv3bqqrq7l+/TqzZ8/m+eefZ+HChdrWR6KWipeIiIwpwWCQ2tpaamtr6e3tZf78+VRUVDB37lwVLol6Kl4iIjImXLt2jerqanbv3k1/fz95eXlUVFSQk6ObKCR2qHiJiIinuru7CQQCNDY2EgqFWLZsGT6fjxkzZngdTSTsVLxERMQTFy9eJBAIsG/fPgBWrlyJz+dj8uTJHicTcY+Kl4iIRNTZs2dxHIeWlhaSkpJYs2YNZWVlTJw40etoIq5T8RIRkYg4deoUjuNw5MgRUlNT8fl8lJSUMH78eK+jiUSMipeIiLjGWsuJEydwHIe2tjbS09N58sknWbt2LWlpaV7HE4k4FS8REQk7ay0nT57k3/7t3+jo6GDChAls2LCBwsJCUlJSvI4n4hkVLxERCZtQKMSBAwfw+/2cP3+eSZMm8eKLL7Jy5UqSkvQrR0T/CkRE5JENDg7S1NREIBDg8uXLTJs2jYqKCtavX69tfUQ+RMVLREQe2sDAAA0NDVRXV3PlyhVmzZrF66+/zuLFi2lvb1fpErmNipeIiIxab28v9fX11NTUcOPGDR577DE2b97M448/rm19RO5BxUtERB7YjRs3qKmpoa6ujr6+PhYsWIDP52Pu3LleRxOJCipeIiJyX1euXKGqqoo9e/YwMDDAkiVLqKioYNasWV5HE4kqKl4iInJXly9fxu/309TURCgUIj8/n/LycqZNm+Z1NJGopOIlIiIfc/78eQKBAM3NzSQkJLBq1SrKy8uZNGmS19FEopqKl4iI3NLZ2YnjOBw6dIjk5GSKi4spKytjwoQJXkcTiQkqXiIiQnt7O47jcOzYMdLS0li3bh3FxcWMGzfO62giMcW14mWM+SKwFTgNvA68CxQCP7bWvu3WeUVE5MFYazl69Ch+v5+TJ08yfvx4nn76adasWUNqaqrX8URikivFyxgzG1hirV1jjPlPwB8CaUAR0GmM+bq1ttONc0t82tnYwfZdrXR2B8nOOs62jXlUFuR4HUtkTLLWcvDgQRzH4ezZs0ycOJHnnnuOgoICkpOTvY4nEtPcuuL1DNBjjPlX4DjQArwP5AN+oABQ8ZKw2NnYwVs7mgkODAHQ0R3krR3NACpfIh8yNDTE/v378fv9XLx4kcmTJ/PSSy+Rn59PYmKi1/FE4oJbxWsmMMdau8EY82VgEtAGTAc+ADJdOq/Eoe27Wm+VrpuCA0Ns39Wq4iXC8D6KjY2NVFVV0d3dzYwZM9i6dStLly7Vlj4iEeZW8boG/GLk418w/BZjKsNXvgqBE7c/wRjzJvAmQE5ODm1tbS5Fu7Ourq6Ink/Cp7M7eNfjkf4+kvDQv8fwGBgY4PDhwxw4cIBgMMi0adN46qmnmD17NsYYTp486er5NY6xQeMYXm4Vr2qG53UBrAb6gWJr7feNMeuAb9/+BGvte8B7AEVFRTY3N9elaHfnxTnl0WVnHafjDuUrOytdYxrFNHYPLxgMUldXR21tLcFgkHnz5lFRUUFubm7E91HUOMYGjWP4uFK8rLUNxphOY0w1cBh4A/iaMaYW2KWJ9RJO2zbmfWSOF0B6ciLbNuZ5mEok8q5du0ZNTQ319fX09/ezaNEiKioqmD17ttfRRGSEa8tJWGu/cNuhN9w6l8S3m/O4fnlXY7ruapS40tPTQyAQoLGxkaGhIZYtW4bP52PGjBleRxOR22gBVYkJlQU5VBYMzw3UJXGJF11dXfj9fvbt2wdAfn4+Pp+PKVOmeJxMRO5GxUtEJMqcO3cOx3FoaWkhMTGRoqIiysrKyMzUDeMiY52Kl4hIlDh9+jSO43D48GFSUlIoKyujpKSEjIwMr6OJyAMaVfEyxiRaa4fu/0gREQkHay0nTpzA7/dz4sQJ0tPTWb9+PWvXriU9Pd3reCIySqO94vUfjTGXrLX/tytpREQEGC5chw8fxnEcOjo6yMjIYMOGDRQWFpKSkuJ1PBF5SA9cvIwx04GXgAr34oiIxLdQKERLSwuO43D+/HmysrJ44YUXWLVqFUlJmh0iEu1G86/4beAPrLW9boUREYlXQ0NDNDU1EQgEuHTpElOnTqWyspIVK1ZoWx+RGHLP4mWM+eTIh+MY3vi6zhgz70MP+Z61ts+tcCIisW5gYIA9e/ZQVVXFlStXmDlzJq+99hpLliyJ+CrzIuK++13xmvChj7962+cA+qkgIvIQ+vr6qK+vp6amhuvXr/PYY4+xefNmHn/8cRUukRh2z+KlSfQiIuF148YNampqqKuro6+vj8cff5yKigrmzp3rdTQRiQDN1BQRiYCrV69SVVVFQ0MDAwMDLF68mIqKCrKzs72OJiIRpOIlIuKiy5cvEwgE2Lt3L6FQiBUrVlBeXs706dO9jiYiHlDxEhFxwYULF/D7/TQ3N5OQkMCqVasoLy9n0qRJXkcTEQ89VPEyxuQAv2utfSvMeUREotqZM2dwHIeDBw+SnJxMcXExpaWlTJw40dXz7mzsYPuuVjq7g2RnpbNtYx6VBTlhe7yIhMd9i5cxZiXQYq0dGPk8C/gx8CcuZxMRiRrt7e04jsOxY8dITU2loqKCkpISxo0b5/q5dzZ28NaOZoIDwzu6dXQHeWtHM8Ady9RoHy8i4XO/dbwM8BPggjHmb4F3gX8A/sxa+6MI5BMRGbOstRw7dgzHcTh58iTjxo3j6aefpqioiLS0tIjl2L6r9VaJuik4MMT2Xa13LFKjfbyIhM/9lpOwxpg24Gng94DDwHettT+MQDYRkTHJWsuhQ4dwHIczZ84wceJENm3axOrVq0lOTo54ns7uoKvHRSR8HmSOl7HWBoEvG2M+AP7aGPOutfasu9FERMaWUChEc3Mzfr+fixcvMnnyZDZv3szKlStJTEz0LFd2VjoddyhN2VnpYXm8iITPAxWvmx9Ya+uMMf878CNjzBMjhUxEJKYNDg6yd+9eAoEA3d3dTJ8+na1bt7J06dIxsY/ito15H5mzBZCenMi2jXlhebyIhM+DFK+P/FSx1jYYY77F8KbZuqtRRGJWf38/DQ0NVFVVce3aNXJycti0aROLFi0aU9v63JyX9aB3KY728SISPg9SvJ6+w7G/YXiSvYhIzOnt7aW2tpba2lqCwSC5ubm8/PLLzJs3b0wVrg+rLMgZVXEa7eNFJDzuW7ystVfvcGwI+IQriUREPHLt2jVqamqor6+nv7+fRYsW4fP5mDNnjtfRRCRGaOV6EYl7PT09VFVVsWfPHgYHB1m2bBk+n4+ZM2d6HU1EYswDFy9jzPvW2ifdDCMiEkldXV34/X727dsHQH5+PuXl5UydOtXjZCISqx5k5fqvWGu/AERuNUARERedO3cOv9/PgQMHSEhIoLCwkLKyMrKysryOJiIx7kGueJW6nkJEJAJOnz6N4zgcPnyYlJQUSktLKS0tJSMjw+toIhInRjPHK9UY82fAfmCntbbfpUwiIvf1oJs8W2tpa2vDcRxOnDhBWloaTzzxBMXFxaSna8FQEYms0RSvQaAJyAe+aIz5lLW22Z1YIiJ39yCbPFtrOXLkCI7jcPr0aTIyMnj22WcpLCwkNTXVs+wiEt9GU7yGrLU/AH5gjPkH4J+NMWXW2l6XsomI3NG9Nnl+aeUsDh48iOM4nDt3jszMTJ5//nkKCgpIStKN3CLirQf5KXRk5E9784C19qgx5u+A5cBuN4KJiNzNnTZzNoQYd/UU7777Ll1dXUyZMoUtW7awYsUKT/dRFBH5sAdZQPVTIx+a246/60oiEZH7+PAmz4mEWJR4keVJZ8lI6Cc5eSavvfYaixcvHhP7KIqIfNhorru/7FoKEZFQsiiyAAAeCklEQVRR2LYxj/+0o5G59izLks4xzgxywWawovhJfmNj8Zjd1kdE5EHW8cq31u6z1p697XgRsNdaO+haOhGR29y4cYPM7sP8SnozQwP9dAxN5EDaXD6zaS0vr57tdTwRkXt6kCte/xV4yhjzZ9batwGMMXOB7wPrgFO3P8EYswz4CXBy5NBvAv8BKAR+fPN1REQe1NWrV6murmb37t0MDAywePFifD4fOTna6FlEosdo3mosBTDGTAC+B3zGWvux0jUiA/imtfZPR55TyvDK90VApzHm69bazodOLSJxo7u7m0AgQGNjI6FQiOXLl+Pz+Zg+fbrX0URERm00xcsaY2YA/wD8R2vtz+/x2AzgOWPMs8AF4P2R//IBP1AAqHiJyF1dvHjx1j6KxhhWrVpFeXk5kydP9jqaiMhDu2fxMsZsBCYbYzYAUxkuT/8nkDRyDGvtv97hqReBb1lrv2aMeQfIAk4A04EPgMyw/Q1EJKacOXOGDz74gPb2dpKSkli7di1lZWVMnDjR62giIo/sfle8VjF89apg5M/JwEKGixQMr+31seJlrW1ieJV7gFrg74HPAS0Mz/M6cftzjDFvAm8C5OTk0NbWNrq/ySPq6uqK6PnEHRrH6HX+/Hn27dtHR0cHSUlJrFixgqVLl5KWlsalS5e4dOmS1xFllPTvMTZoHMPrnsXLWvsXxpiNI38+AzwPfBX4G2vtzrs9zxjza0CmtfbrQAnwFaDYWvt9Y8w64Nt3ONd7wHsARUVFNjc392H/Tg/Ni3NK+Gkco4e1luPHj+M4Du3t7YwbN46nnnqKGTNmsGjRIq/jSRjo32Ns0DiGz6j2zxhZsb4S+CdjTJa19lt3eegPGd5S6NNAO/DrwNeMMbXALk2sF4lv1lpaW1txHIfOzk4mTJjAxo0bWb16NSkpKRG/4i0iEikPUrxurkTYCmCt7TPGfAJ43xhz1Frrv/0JI/s3vnTb4TceKamIRL1QKMT+/fvx+/1cuHCBSZMm8eKLL7Jy5UrtoygiceFBftL9ZwBr7e/ePGCt7TXGfBLQBtkicl+Dg4M0NTURCAS4fPky06ZN45VXXmHZsmXa1kdE4sr97mo8y/CVre8xPIl+FbCXX14FewwodzWhiESt/v5+9uzZQ1VVFVevXiU7O5sNGzaQl5enbX1EJC7d74pXi7X214wx9dbaTxhjtltrt938n8aYKpfziUgU6u3tpa6ujtraWm7cuEFubi5btmxh/vz5KlwiEtfuV7xyRpZ5mDby58qRtxi7gH3AJ9wOKCLR4/r169TU1FBfX09fXx8LFy6koqKCOXPmeB1NRGRMuF/x+j1gPvAphtfs+nPgMwy/7fgfgHnGmFdHJtOLSJy6cuUKVVVVNDQ0MDg4yNKlS/H5fMyaNcvraCIiY8r9itdrwCJ+uRiqAZ4E3gX+PfCbKl0i8evSpUv4/X6ampqw1pKfn4/P52Pq1KleRxMRGZPuV7y+zHD5eg74AsMr2BcAk4DXrbVn3I0XHXY2drB9Vyud3UGys9LZtjGPyoIcr2OJuOb8+fP4/X72799PQkICq1evpry8nKysrPs/WUQkjt2veC1geH/FcQy/5ZgNpI98XgDEffHa2djBWzuaCQ4MAdDRHeStHc0AKl8Sczo6OnAch9bWVpKTkykpKaG0tJQJEyZ4HU1EJCrcr3itAHKAbwJTgPMM79n4AfBfRiba/5q11roZcizbvqv1Vum6KTgwxPZdrSpeEhOstbS3t+M4DsePHyctLY1169ZRXFzMuHHjvI4nIhJV7le8/hqYY609fvOAMeaatbYL+PfGmLnxXLoAOruDozouEi2stRw5cgS/38+pU6cYP348zzzzDEVFRaSmpnodT0QkKt1vk+wB4Phtx779oY/bXcoVNbKz0um4Q8nKzkr3II3IowuFQhw8eBC/38/Zs2fJzMzkueeeo6CggOTkZK/jiYhENW2O9oi2bcz7yBwvgPTkRLZtzPMwlcjoDQ0N0dzcjN/vp6uriylTprBlyxZWrFhBYmKi1/FERGKCitcjujmPS3c1SrQaGBhg7969BAIBenp6mDFjBq+++ipLlizRPooiImGm4hUGlQU5KloSdfr6+ti9ezfV1dVcv36d2bNn8/zzz7Nw4UJt6yMi4hIVL5E4EwwGqa2tpba2lt7eXubPn09FRQVz585V4RIRcZmKl0icuHbtGtXV1ezevZv+/n7y8vKoqKggJ0dXa0VEIkXFSyTGdXd3EwgEaGxsJBQKsXz5csrLy5kxY4bX0URE4o6Kl0iMunjxIn6/n+bm4Z0UVq5cic/nY/LkyR4nExGJXypeIjHm7NmzOI5DS0sLSUlJrFmzhtLSUjIzM72OJiIS91S8RGLEqVOncByHI0eOkJqais/no6SkhPHjx3sdTURERqh4iUQxay3Hjx/HcRza29tJT0/nySefZO3ataSlpXkdT0REbqPiJRKFrLW0trbiOA6dnZ1MmDCBDRs2UFhYSEpKitfxRETkLlS8RKJIKBTiwIED+P1+zp8/T1ZWFi+++CIrV64kKUn/nEVExjr9pBaJAkNDQzQ1NeH3+7l8+TLTpk3j5ZdfZvny5drWR0Qkiqh4iYxhAwMDNDQ0UF1dzZUrV5g1axavv/46ixcv1irzIiJRSMVLZAzq7e2lvr6empoabty4wdy5c3nppZeYP3++CpeISBRT8RIZQ27cuEFNTQ11dXX09fWxYMECKioqeOyxx7yOJiIiYaDiJTIGXLlyherqahoaGhgYGGDp0qX4fD5mzZrldTQREQkjFS8RD126dIlAIEBTUxOhUIj8/HzKy8uZNm2a19FERMQFKl4iHjh//jx+v5/9+/eTkJBAQUEBZWVlTJo0yetoIiLiIhUvkQjq7OzEcRwOHTpEcnIyJSUllJaWMmHCBK+jiYhIBKh4ibjMWkt7ezt+v59jx46RlpbGunXrKC4uZty4cV7HExGRCFLxEnGJtZajR4/iOA6nTp1i/PjxPP3006xZs4bU1FSv44mIiAdUvETCzFrLwYMHcRyHs2fPMnHiRJ577jkKCgpITk72Op6IiHjI1eJljHkV+CqQDbwHFAI/tta+7eZ5RbwwNDTE/v378fv9XLx4kcmTJ/PSSy+Rn59PYmKi1/FERGQMcK14GWNmAr8GnAJKgDSgCOg0xnzdWtvp1rnj1c7GDrbvaqWzO0h2VjrbNuZRWZDjdayYNzg4yN69ewkEAnR3dzNjxgy2bt3K0qVLtY+ii/T9LiLRyM0rXl8B/gD4HsNXut4H8gE/UACoeIXRzsYO3trRTHBgCICO7iBv7WgG0C8jl/T397N7926qq6u5du0aOTk5bNq0iUWLFmlbH5fp+11EopUrxcsY81ngX6217SO/gDKBk8B04IORzyWMtu9qvfVL6KbgwBDbd7XqF1GYBYNB6urqqK2tJRgMMm/ePF555RVyc3NVuCJE3+8iEq3cuuK1BcgwxvwmsBRYCPw20MLw1a8Ttz/BGPMm8CZATk4ObW1tLkW7s66uroieL9w6u4N3PR7pr6WX3BzHYDBIS0sLhw4dYnBwkNmzZ5Ofn39rlfn29nbXzh1v7jeO+n6PDtH+c1WGaRzDy5XiZa198ebHxpga4P8Atlprv2+MWQd8+w7PeY/hCfgUFRXZ3NxcN6LdkxfnDJfsrON03OGXUXZWelT/vR5GuP++PT09BAIBGhsbGRoaYtmyZfh8PmbMmBHW88jt87Z67jpvS9/v0UPjERs0juETqeUkaoDPGGNqgV2aWB9+2zbmfWTOC0B6ciLbNuZ5mCq6dXV14ff72bdvHwD5+fn4fD6mTJnicbLYNJp5W/p+F5Fo5XrxstaWjHz4htvnimc3fzHpLq9Hd+7cORzHoaWlhcTERIqKiigrKyMzU1MT3TSaeVv6fheRaKUFVGNIZUGOfvE8gtOnT+M4DocPHyYlJYWysjJKSkrIyMjwOlpcuNe8rTvR97uIRCMVL4lr1lpOnDiB4zi0tbWRnp7O+vXrWbt2Lenp6V7HiyvZWel3nbclIhIrVLwkLllrOXz4MI7j0NHRQUZGBhs2bKCwsJCUlBSv48UlzdsSkXig4iVxJRQK0dLSguM4nD9/nqysLF544QVWrVpFUpL+OXhJ87ZEJB7oN43EhaGhIZqamggEAly6dImpU6dSWVnJihUrtK3PGHJz3lZbW5tuXxeRmKTiJTFtYGCAPXv2UFVVxZUrV5g1axavv/46ixcv1irz96G9EEVEwk/FS2JSX18f9fX11NTUcP36dR577DE2b97M448/rsL1ALQXooiIO1S8JKb09vby85//nLq6Ovr6+liwYAE+n4+5c+d6HS2qaC9EERF3qHhJTLh69SpVVVXs3r2bwcFBlixZgs/nIzs72+toUWm0a2qJiMiDUfGSqHb58mUCgQB79+4lFAoxb948Nm3adGvjank4WlNLRMQdKl4SlS5cuIDf76e5uZmEhARWrVpFeXk5PT09Kl1hoDW1RETcoeIlUaWzsxO/38/BgwdJTk6muLiY0tJSJk6cCEBPT4/HCWOD1tQSEXGHipdEhfb2dhzH4dixY6SlpbFu3TqKi4sZN26c19FilvZCFBEJPxUvGbOstRw7dgzHcTh58iTjxo3j6aefZs2aNaSmpnodT0REZNRUvCTsHnXhTWsthw4dwnEczpw5w8SJE9m0aROrV68mOTnZxeRjmxY0FRGJfipeElaPsvBmKBSiubkZv9/PxYsXmTx5Mps3b2blypUkJia6nn0s04KmIiKxQcVLwuphFt4cHBxk7969BAIBuru7mT59Olu3bmXp0qXaR3GEFjQVEYkNKl4SVqNZeLO/v5+Ghgaqq6u5evUqOTk5bNq0iUWLFmlbn9toQVMRkdig4iVh9SALb/b29lJXV0dNTQ3BYJDc3FwqKyuZN2+eCtddaEFTEZHYoOIlYXWvhTevX79OdXU19fX19Pf3s3DhQioqKpgzZ46HiaODFjQVEYkNKl4SVndaePPz63JIPbuP//r/7WFwcJBly5bh8/mYOXOmx2mjhxY0FRGJDSpeEnY3F968dOkSfr+fpv/5zwDk5+dTXl7O1KlTPU4YnbSgqYhI9FPxkrA7d+4cfr+fAwcOkJCQQGFhIWVlZWRlZXkdTURExFMqXhI2p0+fxnEcDh8+TEpKCqWlpZSWlpKRkeF1NBERkTFBxUseibWWtrY2HMfhxIkTpKens379etauXUt6uu64ExER+TAVL3ko1lqOHDmC4zicPn2ajIwMnn32WYqKikhJSfE6noiIyJik4iWjEgqFaGlpwe/3c+7cOTIzM3n++ecpKCggKUnfTiIiIvei35TyQIaGhti3bx9+v59Lly4xdepUKisrWb58edzvoygiIvKgVLzkngYGBmhsbCQQCHDlyhVmzpzJa6+9xuLFi7WPooiIyCipeMkd9fX1UV9fT01NDdevX2fOnDm8+OKLLFiwQNv6iIiIPCQVL/mIGzduUFtbS11dHb29vcyfP5+Kigrmzp2rwiUiIvKIVLwEgKtXr1JdXc3u3bsZGBhg8eLF+Hw+cnK0UrqIiEi4qHjFue7ubgKBAI2NjYRCIZYvX47P52P69OleRxMREYk5Kl5x6uLFi/j9fvbt24cxhlWrVlFeXs7kyZO9jiYiIhKzXClexpgM4J+BiUAf8O+A3wFeAPYCb1hrrRvnlns7c+YMfr+flpYWkpKSWLt2LWVlZUycONHraCIiIjHPrSteWcAXrbX7jDF/BGwEnrDWFhpjdgPFQI1L55Y7OHnyJI7jcPToUVJTU6moqKC4uJjx48d7HU1ERCRuuFK8rLWnjTHjjTENwFXgK8AvjDEzgFagEBUv11lrOX78OI7j0N7ezrhx43jqqadYs2YNaWlpXscTERGJO67N8bLWtgKFxpg/ASYDF4GVwH8HVrh1XhkuXK2trTiOQ2dnJxMmTGDjxo0UFhaSnJzsdTwREZG45dYcr3nANWvtBeBHwA6Gr3qdBZYBPXd4zpvAmwA5OTm0tbW5Ee2uurq6Ino+N4RCIdra2mhubqa7u5uMjAxKS0t5/PHHSUxMpKOjw+uIrouFcRSNY6zQOMYGjWN4uXXFqxhYAvwJsBr4b4APeJfhSfbfuv0J1tr3gPcAioqKbG5urkvR7s6Lc4bD4OAgTU1NBAIBLl++zLRp03jllVdYtmxZ3Gzrs7Oxg+27WunsDpKd1cO2jXlUFmgNsmgWrf8e5aM0jrFB4xg+bhWvHwDfMMb8L6AL+BTwhwzP62oCal06b1zp7++noaGB6upqrl69SnZ2Nhs2bCAvLy+uVpnf2djBWzuaCQ4MAdDRHeStHc0AKl8iIjKmuDW5fgD45G2H3x75Tx5Rb28vdXV11NbWcuPGDXJzc6msrGTevHlxVbhu2r6r9Vbpuik4MMT2Xa0qXiIiMqZoAdUocv36dWpqaqivr6evr4+FCxdSUVHBnDlzvI7mqc7u4KiOi4iIeEXFKwpcuXKFQCDAnj17GBwcZOnSpVRUVDBz5kyvo40J2VnpdNyhZGVnpXuQRkRE5O5UvMawS5cu4ff7aWpqwlpLfn4+Pp+PqVOneh1tTNm2Me8jc7wA0pMT2bYxz8NUIiIiH6fiNQadP38ex3E4cOAACQkJrF69mvLycrKysryONibdnMf1y7sa03VXo4iIjEkqXmNIR0cHjuPQ2tpKSkoKpaWllJSUMGHCBK+jjXmVBTlUFgyv/6bbnkVEZKxS8fKYtZb29nYcx+H48eOkpaXxxBNPUFxcTHq65iiJiIjEEhUvj1hrOXr0KI7jcOrUKcaPH88zzzxDUVERqampXscTERERF6h4RVgoFOLgwYP4/X7Onj1LZmYmzz33HAUFBdpHUUREJMapeEXI0NAQzc3N+P1+urq6mDJlClu2bGHFihUkJiZ6HU9EREQiQMXLZQMDA+zdu5dAIEBPTw8zZszg1VdfZcmSJXGzj6KIiIgMU/FySV9fH7t376a6uprr168zZ84cXnjhBRYsWBCX2/qIiIiIilfYBYNBamtrqa2tpbe3l/nz51NRUcHcuXNVuEREROKcileYXLt2jerqanbv3k1/fz95eXlUVFSQk6NFPEVERGSYitcj6u7uJhAI0NjYSCgUYvny5ZSXlzNjxgyvo4mIiMgYo+L1kC5evEggEGDfvn0ArFy5Ep/Px+TJkz1OJiIiImOVitconT17FsdxaGlpISkpiaKiIsrKysjMzPQ6moiIiIxxKl4P6NSpUziOw5EjR0hNTcXn81FSUsL48eO9jiYiIiJRQsXrHqy1nDhxAsdxaGtrIz09nSeffJK1a9eSlpbmdTwRERGJMiped2CtpbW1Fb/fT0dHBxMmTGDDhg0UFhaSkpLidTwRERGJUipeHxIKhThw4AB+v5/z588zadIkXnzxRVauXElSkr5UIiIi8mjUJoDBwUEOHz7Mj370Iy5fvsy0adN4+eWXWb58ubb1ERERkbCJ++J148YN/uZv/oYrV64wa9YsXn/9dRYvXqxV5kVERCTs4r54jRs3juXLlzNu3DjKyspUuMa4nY0dbN/VSmd3kOysdLZtzKOyQLsDiIhIdIj74gXw7LPP0tbWptI1xu1s7OCtHc0EB4YA6OgO8taOZgCVLxERiQqawCRRY/uu1lul66bgwBDbd7V6lEhERGR0VLwkanR2B0d1XEREZKxR8ZKokZ2VPqrjIiIiY42Kl0SNbRvzSE9O/Mix9OREtm3M8yiRiIjI6GhyvUSNmxPodVejiIhEKxUviSqVBTkqWiIiErX0VqOIiIhIhKh4iYiIiESIipeIiIhIhLhSvIwxqcaY7xpj/MaYH418/rfGmD3GmD9z45wiIiIiY51bV7xeA/Zaa33AEeCLQBpQBLxpjMl26bwiIiIiY5ZbxasJ+O7Ix4OABd4H8gE/UODSeUVERETGLFeWk7DWNgMYYyqBTOAKcBGYDnwwckxEREQkrri2jpcx5teBVcDngP8NSAVagELgxB0e/ybwJkBOTg5tbW1uRbujrq6uiJ5P3KFxjA0ax9igcYwNGsfwcqV4GWPmAluttS+PfN4w8vn3jTHrgG/f/hxr7XvAewBFRUU2NzfXjWgfs7Ox40MrofdoJfQYEKnvHXGXxjE2aBxjg8YxfNy64vUZYIUxxj/y+d8CWcaYWmCXtbbTpfOOys7GDt7a0UxwYAiAju4gb+1oBlD5EhERkbBza47X28Dbtx3+2FUur23f1XqrdN0UHBhi+65WFS8REREJu7heQLWzOziq4yIiIiKPIq6LV3ZW+qiOi4iIiDyKuC5e2zbmkZ6c+JFj6cmJbNuY51Gi6LOzsYPyL/2ceV/8fyn/0s/Z2djhdSQREZExy7XlJKLBzXlcv7yrMV13NY6Cbk4QEREZnbguXjBcECoLhtcN0+2yo6ObE0REREYnrt9qlEejmxNERERGR8VLHppuThARERkdFS95aLo5QUREZHTifo6XPDzdnCAiIjI6Kl7ySG7enCAiIiL3p7caRURERCJExUtEREQkQlS8RERERCJExUtEREQkQlS8RERERCJExUtEREQkQlS8RERERCJExUtEREQkQlS8RERERCJExUtEREQkQlS8RERERCJExUtEREQkQoy11usMH2OMuQC0R/i0U4GLET6nhJ/GMTZoHGODxjE2aBzvb661dtqDPHBMFi8vGGN2W2uLvM4hj0bjGBs0jrFB4xgbNI7hpbcaRURERCJExUtEREQkQlS8fuk9rwNIWGgcY4PGMTZoHGODxjGMNMdLREREJEJ0xUtEREQkQuK+eJlhf2uM2WOM+TOv88jDMcakGmO+a4zxG2N+ZIxJ9TqTPBxjzKvGmDNe55BHY4z5ojGm3hjzQ2NMktd5ZHRGfqb+48jP1J8YY2Z4nSlWxH3xAkqANKAIeNMYk+1xHnk4rwF7rbU+4Aiw0eM88hCMMTOBXwNOeZ1FHp4xZjawxFq7BtgDPOtxJBm9Z4BTIz9Ta4EKj/PEDBUvKATeB/IBP1DgbRx5SE3Ad0c+HkSL/UWrrwB/AIS8DiKP5Bmgxxjzr0AOsMvjPDJ6bcDvG2NagU3ATm/jxA4VL8hk+Jf0dOCDkc8lylhrm621p40xlUCmtbbK60wyOsaYzwL/aq2N9K4VEn4zgTnW2g3AFWCrx3lk9M4AL1tr84D/B/i0t3Fih4oX9ACpwHFgwsjnEoWMMb8O+IDPeZ1FHsoW4DeNMR8AS40xX/E4jzy8a8AvRj7+BbDUwyzycLYxXL5g+GpXuYdZYoqKFzQAxdbao8A6oNHjPPIQjDFzga3W2j+0WiMlKllrX7TWrrfWrgdarLVf8DqTPLRqoHjk49VAq4dZ5OGVjfxZChzzMkgsUfGCGiDLGFML1FtrO70OJA/lM8CKkTtw/MaYT3odSCReWWsbgE5jTDWwAPgfHkeS0ftr4JWRMdwM/F8e54kZWkBVREREJEJ0xUtEREQkQlS8RERERCJExUtEREQkQlS8RERERCJExUtEREQkQlS8RCSqGGN+2xjzRphf864/C40xv2uM+XQ4zyci8Us7xotIVDDGpFlrez/0+WQg2Vp77kPH/h7+//buJ8SmMI7D+POLDUVslbIQK1IiK1Ps1CBhLYUFipXUbGyGpbJjo5SyYUo2SpKFLNXESkKWhJrRpOZrMefW6bpzx0zjmsnzqbt4//a7m9O3955zLlu6lr5IcrE1Zy3wBEjzmQY2VdWBJL5AWdJfZfCStFw8rqrhVnsjMAIc63QkmfPFuUm+A7s67apaDbwEXvdZtmre1UpSDwYvSUteVW0AJpJ8qyoAkryqqq9VtQ94AzxoLVnBzPVtqmlPJRmaZfvzwK0kU7OMA1yoqm3AlfYJmyTNl8FL0nJwGbjdo/8acDPJfmBPp7OqjgI7koz027SqdgNn+P3nyW5XgXfAo6p6CIwm+fnn5UvSDG+ul7SkNX+AvjXJve6xJG+Bz82JWNt2YHyOfXcC14G7wFhVre83P8kzYC/wxdAlaaE88ZK0pCV5X1WHW13VNX683a6qzcBBYLTXflW1EjgNnAQOJflUVUeAp1U1nORjn1omgRsL+yaSZPCStAwkmayqNcBzYB1wqte8qjrbjJ1oPwHZ5Q7wARhKMtHsf7+qpoFzwKXFrl+SOirJv65BkhZFc8/WeHMyNducyjwufM1Tj0nyYzFqlPR/M3hJkiQNiDfXS5IkDYjBS5IkaUAMXpIkSQNi8JIkSRoQg5ckSdKAGLwkSZIG5Bf/Zg06/yEFQAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "poly_fit = np.polyfit(x, y, 1)\n",
    "poly_1d = np.poly1d(poly_fit)\n",
    "xs = np.linspace(x.min(), x.max())\n",
    "ys = poly_1d(xs)\n",
    "\n",
    "fig = plt.figure(figsize=(10, 6))\n",
    "ax = fig.add_subplot(111)\n",
    "ax.set_xlabel('小テスト')\n",
    "ax.set_ylabel('期末テスト')\n",
    "ax.plot(xs, ys, color='gray', \n",
    "        label=f'{poly_fit[1]:.2f}+{poly_fit[0]:.2f}x')\n",
    "ax.scatter(x, y)\n",
    "ax.legend()\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 回帰分析における仮説"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### statsmodelsによる回帰分析"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.253141Z",
     "start_time": "2018-08-24T07:49:06.234941Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.676</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.658</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   37.61</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Fri, 24 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>8.59e-06</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>16:49:06</td>     <th>  Log-Likelihood:    </th> <td> -76.325</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   156.7</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    18</td>      <th>  BIC:               </th> <td>   158.6</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     1</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   23.6995</td> <td>    4.714</td> <td>    5.028</td> <td> 0.000</td> <td>   13.796</td> <td>   33.603</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>      <td>    6.5537</td> <td>    1.069</td> <td>    6.133</td> <td> 0.000</td> <td>    4.309</td> <td>    8.799</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 2.139</td> <th>  Durbin-Watson:     </th> <td>   1.478</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.343</td> <th>  Jarque-Bera (JB):  </th> <td>   1.773</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.670</td> <th>  Prob(JB):          </th> <td>   0.412</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.422</td> <th>  Cond. No.          </th> <td>    8.32</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.676\n",
       "Model:                            OLS   Adj. R-squared:                  0.658\n",
       "Method:                 Least Squares   F-statistic:                     37.61\n",
       "Date:                Fri, 24 Aug 2018   Prob (F-statistic):           8.59e-06\n",
       "Time:                        16:49:06   Log-Likelihood:                -76.325\n",
       "No. Observations:                  20   AIC:                             156.7\n",
       "Df Residuals:                      18   BIC:                             158.6\n",
       "Df Model:                           1                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     23.6995      4.714      5.028      0.000      13.796      33.603\n",
       "小テスト           6.5537      1.069      6.133      0.000       4.309       8.799\n",
       "==============================================================================\n",
       "Omnibus:                        2.139   Durbin-Watson:                   1.478\n",
       "Prob(Omnibus):                  0.343   Jarque-Bera (JB):                1.773\n",
       "Skew:                           0.670   Prob(JB):                        0.412\n",
       "Kurtosis:                       2.422   Cond. No.                         8.32\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "formula = '期末テスト ~ 小テスト'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 回帰係数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.257552Z",
     "start_time": "2018-08-24T07:49:06.254216Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 1. ,  4.2],\n",
       "       [ 1. ,  7.2],\n",
       "       [ 1. ,  0. ],\n",
       "       [ 1. ,  3. ],\n",
       "       [ 1. ,  1.5],\n",
       "       [ 1. ,  0.9],\n",
       "       [ 1. ,  1.9],\n",
       "       [ 1. ,  3.5],\n",
       "       [ 1. ,  4. ],\n",
       "       [ 1. ,  5.4],\n",
       "       [ 1. ,  4.2],\n",
       "       [ 1. ,  6.9],\n",
       "       [ 1. ,  2. ],\n",
       "       [ 1. ,  8.8],\n",
       "       [ 1. ,  0.3],\n",
       "       [ 1. ,  6.7],\n",
       "       [ 1. ,  4.2],\n",
       "       [ 1. ,  5.6],\n",
       "       [ 1. ,  1.4],\n",
       "       [ 1. ,  2. ]])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X = np.array([np.ones_like(x), x]).T\n",
    "X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.262050Z",
     "start_time": "2018-08-24T07:49:06.258677Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(23.699, 6.554)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "beta0_hat, beta1_hat = np.linalg.lstsq(X, y)[0]\n",
    "beta0_hat, beta1_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.265451Z",
     "start_time": "2018-08-24T07:49:06.263014Z"
    }
   },
   "outputs": [],
   "source": [
    "y_hat = beta0_hat + beta1_hat * x\n",
    "eps_hat = y - y_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.270016Z",
     "start_time": "2018-08-24T07:49:06.266359Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "134.290"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s_var = np.var(eps_hat, ddof=p+1)\n",
    "s_var"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.274116Z",
     "start_time": "2018-08-24T07:49:06.271523Z"
    }
   },
   "outputs": [],
   "source": [
    "C0, C1 = np.diag(np.linalg.pinv(np.dot(X.T, X)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.278932Z",
     "start_time": "2018-08-24T07:49:06.275370Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4.714, 1.069)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.sqrt(s_var * C0), np.sqrt(s_var * C1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.285336Z",
     "start_time": "2018-08-24T07:49:06.280108Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(13.796, 33.603)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rv = stats.t(n-2)\n",
    "\n",
    "lcl = beta0_hat - rv.isf(0.025) * np.sqrt(s_var * C0)\n",
    "hcl = beta0_hat - rv.isf(0.975) * np.sqrt(s_var * C0)\n",
    "lcl, hcl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.291444Z",
     "start_time": "2018-08-24T07:49:06.286385Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4.309, 8.799)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rv = stats.t(n-2)\n",
    "\n",
    "lcl = beta1_hat - rv.isf(0.025) * np.sqrt(s_var * C1)\n",
    "hcl = beta1_hat - rv.isf(0.975) * np.sqrt(s_var * C1)\n",
    "lcl, hcl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.295776Z",
     "start_time": "2018-08-24T07:49:06.292661Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "6.133"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t = beta1_hat / np.sqrt(s_var * C1)\n",
    "t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.300759Z",
     "start_time": "2018-08-24T07:49:06.296921Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.000"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(1 - rv.cdf(t)) * 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.305075Z",
     "start_time": "2018-08-24T07:49:06.301957Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.028"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t = beta0_hat / np.sqrt(s_var * C0)\n",
    "t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.309479Z",
     "start_time": "2018-08-24T07:49:06.306217Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.000"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(1 - rv.cdf(t)) * 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 重回帰モデル"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-24T07:49:06.324298Z",
     "start_time": "2018-08-24T07:49:06.310535Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.756</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.727</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   26.35</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Fri, 24 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>6.19e-06</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>16:49:06</td>     <th>  Log-Likelihood:    </th> <td> -73.497</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   153.0</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    17</td>      <th>  BIC:               </th> <td>   156.0</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     2</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   -1.8709</td> <td>   11.635</td> <td>   -0.161</td> <td> 0.874</td> <td>  -26.420</td> <td>   22.678</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>      <td>    6.4289</td> <td>    0.956</td> <td>    6.725</td> <td> 0.000</td> <td>    4.412</td> <td>    8.446</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>睡眠時間</th>      <td>    4.1917</td> <td>    1.778</td> <td>    2.357</td> <td> 0.031</td> <td>    0.440</td> <td>    7.943</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 2.073</td> <th>  Durbin-Watson:     </th> <td>   1.508</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.355</td> <th>  Jarque-Bera (JB):  </th> <td>   1.716</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.660</td> <th>  Prob(JB):          </th> <td>   0.424</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.437</td> <th>  Cond. No.          </th> <td>    38.0</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.756\n",
       "Model:                            OLS   Adj. R-squared:                  0.727\n",
       "Method:                 Least Squares   F-statistic:                     26.35\n",
       "Date:                Fri, 24 Aug 2018   Prob (F-statistic):           6.19e-06\n",
       "Time:                        16:49:06   Log-Likelihood:                -73.497\n",
       "No. Observations:                  20   AIC:                             153.0\n",
       "Df Residuals:                      17   BIC:                             156.0\n",
       "Df Model:                           2                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     -1.8709     11.635     -0.161      0.874     -26.420      22.678\n",
       "小テスト           6.4289      0.956      6.725      0.000       4.412       8.446\n",
       "睡眠時間           4.1917      1.778      2.357      0.031       0.440       7.943\n",
       "==============================================================================\n",
       "Omnibus:                        2.073   Durbin-Watson:                   1.508\n",
       "Prob(Omnibus):                  0.355   Jarque-Bera (JB):                1.716\n",
       "Skew:                           0.660   Prob(JB):                        0.424\n",
       "Kurtosis:                       2.437   Cond. No.                         38.0\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "formula = '期末テスト ~ 小テスト + 睡眠時間'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 回帰係数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.081894Z",
     "start_time": "2018-08-20T18:22:38.079369Z"
    }
   },
   "outputs": [],
   "source": [
    "x1 = df['小テスト']\n",
    "x2 = df['睡眠時間']\n",
    "y = df['期末テスト']\n",
    "p = 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.088678Z",
     "start_time": "2018-08-20T18:22:38.083316Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-1.871, 6.429, 4.192)"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X = np.array([np.ones_like(x1), x1, x2]).T\n",
    "beta0_hat, beta1_hat, beta2_hat = np.linalg.lstsq(X, y)[0]\n",
    "beta0_hat, beta1_hat, beta2_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.104166Z",
     "start_time": "2018-08-20T18:22:38.089732Z"
    }
   },
   "outputs": [],
   "source": [
    "y_hat = beta0_hat + beta1_hat * x1 + beta2_hat * x2\n",
    "eps_hat = y - y_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.108484Z",
     "start_time": "2018-08-20T18:22:38.105875Z"
    }
   },
   "outputs": [],
   "source": [
    "s_var = np.sum(eps_hat ** 2) / (n - p - 1)\n",
    "C0, C1, C2 = np.diag(np.linalg.pinv(np.dot(X.T, X)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.114358Z",
     "start_time": "2018-08-20T18:22:38.109414Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.440, 7.943)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rv = stats.t(n-p-1)\n",
    "\n",
    "lcl = beta2_hat - rv.isf(0.025) * np.sqrt(s_var * C2)\n",
    "hcl = beta2_hat - rv.isf(0.975) * np.sqrt(s_var * C2)\n",
    "lcl, hcl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ダミー変数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.130955Z",
     "start_time": "2018-08-20T18:22:38.115424Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.782</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.724</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   13.46</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Tue, 21 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>7.47e-05</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>03:22:38</td>     <th>  Log-Likelihood:    </th> <td> -72.368</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   154.7</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    15</td>      <th>  BIC:               </th> <td>   159.7</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     4</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "       <td></td>          <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th>   <td>   -0.4788</td> <td>   12.068</td> <td>   -0.040</td> <td> 0.969</td> <td>  -26.202</td> <td>   25.244</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>通学方法[T.徒歩]</th>  <td>   -5.8437</td> <td>    5.447</td> <td>   -1.073</td> <td> 0.300</td> <td>  -17.453</td> <td>    5.766</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>通学方法[T.自転車]</th> <td>    1.8118</td> <td>    6.324</td> <td>    0.286</td> <td> 0.778</td> <td>  -11.668</td> <td>   15.292</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>        <td>    6.0029</td> <td>    1.033</td> <td>    5.809</td> <td> 0.000</td> <td>    3.800</td> <td>    8.206</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>睡眠時間</th>        <td>    4.5238</td> <td>    1.809</td> <td>    2.501</td> <td> 0.024</td> <td>    0.668</td> <td>    8.380</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 1.764</td> <th>  Durbin-Watson:     </th> <td>   1.418</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.414</td> <th>  Jarque-Bera (JB):  </th> <td>   0.989</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.545</td> <th>  Prob(JB):          </th> <td>   0.610</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.985</td> <th>  Cond. No.          </th> <td>    39.8</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.782\n",
       "Model:                            OLS   Adj. R-squared:                  0.724\n",
       "Method:                 Least Squares   F-statistic:                     13.46\n",
       "Date:                Tue, 21 Aug 2018   Prob (F-statistic):           7.47e-05\n",
       "Time:                        03:22:38   Log-Likelihood:                -72.368\n",
       "No. Observations:                  20   AIC:                             154.7\n",
       "Df Residuals:                      15   BIC:                             159.7\n",
       "Df Model:                           4                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "===============================================================================\n",
       "                  coef    std err          t      P>|t|      [0.025      0.975]\n",
       "-------------------------------------------------------------------------------\n",
       "Intercept      -0.4788     12.068     -0.040      0.969     -26.202      25.244\n",
       "通学方法[T.徒歩]     -5.8437      5.447     -1.073      0.300     -17.453       5.766\n",
       "通学方法[T.自転車]     1.8118      6.324      0.286      0.778     -11.668      15.292\n",
       "小テスト            6.0029      1.033      5.809      0.000       3.800       8.206\n",
       "睡眠時間            4.5238      1.809      2.501      0.024       0.668       8.380\n",
       "==============================================================================\n",
       "Omnibus:                        1.764   Durbin-Watson:                   1.418\n",
       "Prob(Omnibus):                  0.414   Jarque-Bera (JB):                0.989\n",
       "Skew:                           0.545   Prob(JB):                        0.610\n",
       "Kurtosis:                       2.985   Cond. No.                         39.8\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "formula = '期末テスト ~ 小テスト + 睡眠時間 + 通学方法'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## モデルの選択"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.144601Z",
     "start_time": "2018-08-20T18:22:38.131982Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.676</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.658</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   37.61</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Tue, 21 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>8.59e-06</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>03:22:38</td>     <th>  Log-Likelihood:    </th> <td> -76.325</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   156.7</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    18</td>      <th>  BIC:               </th> <td>   158.6</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     1</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   23.6995</td> <td>    4.714</td> <td>    5.028</td> <td> 0.000</td> <td>   13.796</td> <td>   33.603</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>      <td>    6.5537</td> <td>    1.069</td> <td>    6.133</td> <td> 0.000</td> <td>    4.309</td> <td>    8.799</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 2.139</td> <th>  Durbin-Watson:     </th> <td>   1.478</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.343</td> <th>  Jarque-Bera (JB):  </th> <td>   1.773</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.670</td> <th>  Prob(JB):          </th> <td>   0.412</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.422</td> <th>  Cond. No.          </th> <td>    8.32</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.676\n",
       "Model:                            OLS   Adj. R-squared:                  0.658\n",
       "Method:                 Least Squares   F-statistic:                     37.61\n",
       "Date:                Tue, 21 Aug 2018   Prob (F-statistic):           8.59e-06\n",
       "Time:                        03:22:38   Log-Likelihood:                -76.325\n",
       "No. Observations:                  20   AIC:                             156.7\n",
       "Df Residuals:                      18   BIC:                             158.6\n",
       "Df Model:                           1                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     23.6995      4.714      5.028      0.000      13.796      33.603\n",
       "小テスト           6.5537      1.069      6.133      0.000       4.309       8.799\n",
       "==============================================================================\n",
       "Omnibus:                        2.139   Durbin-Watson:                   1.478\n",
       "Prob(Omnibus):                  0.343   Jarque-Bera (JB):                1.773\n",
       "Skew:                           0.670   Prob(JB):                        0.412\n",
       "Kurtosis:                       2.422   Cond. No.                         8.32\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.array(df['小テスト'])\n",
    "y = np.array(df['期末テスト'])\n",
    "p = 1\n",
    "\n",
    "formula = '期末テスト ~ 小テスト'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.149065Z",
     "start_time": "2018-08-20T18:22:38.145691Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 51.225,  70.886,  23.699,  43.361,  33.53 ,  29.598,  36.152,\n",
       "        46.638,  49.914,  59.09 ,  51.225,  68.92 ,  36.807,  81.372,\n",
       "        25.666,  67.61 ,  51.225,  60.4  ,  32.875,  36.807])"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "y_hat = np.array(result.fittedvalues)\n",
    "y_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.153640Z",
     "start_time": "2018-08-20T18:22:38.150159Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 15.775,   0.114,  -4.699,  -8.361,   1.47 ,  10.402, -13.152,\n",
       "        -9.638, -10.914,  -4.09 , -11.225,   1.08 ,  -7.807,   6.628,\n",
       "        21.334,   9.39 ,   0.775,  -5.4  , -14.875,  23.193])"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "eps_hat = np.array(result.resid)\n",
    "eps_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.158021Z",
     "start_time": "2018-08-20T18:22:38.154580Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2417.228"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.sum(eps_hat ** 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 決定係数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.161722Z",
     "start_time": "2018-08-20T18:22:38.159037Z"
    }
   },
   "outputs": [],
   "source": [
    "total_var = np.sum((y - np.mean(y))**2)\n",
    "exp_var = np.sum((y_hat - np.mean(y))**2)\n",
    "unexp_var = np.sum(eps_hat ** 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.166082Z",
     "start_time": "2018-08-20T18:22:38.162816Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(7468.550, 7468.550)"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "total_var, exp_var + unexp_var"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.170409Z",
     "start_time": "2018-08-20T18:22:38.167023Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.676"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "exp_var / total_var"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.174678Z",
     "start_time": "2018-08-20T18:22:38.171375Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.676"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.corrcoef(x, y)[0, 1] ** 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 自由度調整済み決定係数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.178817Z",
     "start_time": "2018-08-20T18:22:38.175704Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.658"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "1 - (unexp_var / (n - p - 1)) / (total_var / (n - 1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### F検定"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.183230Z",
     "start_time": "2018-08-20T18:22:38.179814Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "37.615"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f = (exp_var / p)  / (unexp_var / (n - p - 1))\n",
    "f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.187620Z",
     "start_time": "2018-08-20T18:22:38.184168Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.000"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rv = stats.f(p, n-p-1)\n",
    "1 - rv.cdf(f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 最大対数尤度とAIC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.193573Z",
     "start_time": "2018-08-20T18:22:38.188600Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.031"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prob = 0.3\n",
    "coin_result = [0, 1, 0, 0, 1]\n",
    "\n",
    "rv = stats.bernoulli(prob)\n",
    "L = np.prod(rv.pmf(coin_result))\n",
    "L"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.329480Z",
     "start_time": "2018-08-20T18:22:38.194588Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAl8AAAFoCAYAAABpKJ6VAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xd0leeB5/Hvc1URAhWEDBLNNIOEQQIJRLEQvWNjxbhk7Hg8Ns6Od+YkZyc7u5OzM04mziTrbM7M7MwkaydZm8Q12BDTe1cBCUS3QYBoAoNQA/Xy7B8SWsAIBEh6b/l9zsmJ79Vzr36Xx0I/P+/7Pq+x1iIiIiIincPldAARERERX6LyJSIiItKJVL5EREREOpHKl4iIiEgnUvkSERER6UQqXyIiIiKdSOVLREREpBOpfImIiIh0IpUvERERkU6k8iUiIiLSifydDnA3UVFRdsCAAR36PWprawkMDOzQ7yH3T/PifjQn7knz4n40J+6pM+YlNze3yFrb817j3Lp8DRgwgJycnA79HgUFBXR0wZP7p3lxP5oT96R5cT+aE/fUGfNijDnTlnE67CgiIiLSiVS+RERERDqRypeIiIhIJ1L5EhEREelEKl8iIiIineie5cs0edcYs88Y8+NWxgQZY5YbY/YbY15vfm6UMSbDGLPLGPN+8/vMM8Ycb35ulzEmtL0/kIiIiIg7a8tWEylAMJAEFBpjfm2tLbxtzNPAYeBl4JQx5j3AD3jaWnvJGLMBeAQIBX5qrX2vfeKLiIi4h+rqaq5cuUJFRQXHjh1zOo7cpr6+/oHnJSAggOjoaLp3794uWdpSvsYAW4GRwC4gEbi9fI0B1tBU1LKAwdbafcaYZ40xPwE+ay5hocB3jDFLgKPW2lfb5VOIiIg4qKysjK+//pqePXsSGRlJSEgIxhinY8lNampqCAoKuu/XWWupqqriwoULAO1SwNpyzlcYUAREA9uaH7c2xgLHboyx1n4CxAPDjDH9gTPAv1hrJwABxpgnHvYDiIiIOK2oqIg+ffoQERGBv7+/ipcXMcYQEhJCbGwsly9fbpf3bMvKVxkQBBylaYXrdCtjooADQDpQZowZDRy01tYaY7YDKc1l7IZsYBCw8+Y3al4VWwIQGxtLQUHBfX2g+3X16tUOfX95MJoX96M5cU+aF/dQWVmJy+WipqaG+vp6p+PIHTzsvLhcLqqqqtqll7SlfOUC6dbaPxpjUoH3WxkTR9PK2EggH/gX4BOaDlkmAruNMX8LHLbWrqbpEOWvbn8ja+07wDsASUlJHX5vR0C3gXBTmhf3ozlxT5oX5x07dozg4OCWxw9yeEs63sPOi7+/f7v8vLWlfGUBf2GMyQbW3+Fke4DlwEfAS8BvrbU1zed6/d/m/99lrd1jjDkNfGyM+RGQaa3NfOhPICLtrrGxkbKyMq5fv05jYyONjY1cvHiR+vp6jDF069aN7t27ExQUpMMrIiL36Z7ly1prgZYT440xk4HNwCBr7ZnmMdXAottedx6YcdtzV4BpDx9bRNpLWVkZp06dorCwkJKSEkpKSigtLaWxsfGerw0ICCAsLIywsDB69+5Nv3796NOnD126dOmE5CLS2f7whz/w5ptvkp+fD8DHH39McnIygwYNatPr33vvPX7xi19w+PDhu4579dVXWbZsGeHh4S3PXbp0iVWrVpGVlcXbb79NREREy9eKior43ve+x09+8pOWx7t37+bJJ58kKSmJn/3sZxhjqKmpYe7cuaxZs4bf/OY3fP755/f7R9Au2rLydbscIIFvXvEoIh6gurqakydPcvr0aU6fPk1xcTHQtBwfGRlJr169GD58OJGRkXTv3h2Xy4XL5eLrr78mJiaGxsZGrl27Rnl5OeXl5Vy7do3i4mJ2797Nrl27AOjZsyd9+/blscceY9CgQfj5+Tn5kUWkg6xZs4a//du/JTMzk5iYmFu+Nnjw4G+cH2WtpbGxEX//b9aP//Jf/gs///nPWx5/73vf480332x5PGnSpJZ/fu211/jFL37R8vjll1++5b1qamr49re/TVZWVstzP/jBD3j99dcBOHjwIAEBAW3+nO3tvsuXtbaCpj29RMRDWGspKChg//79HDt2jPr6egIDAxkwYADJyckMHDiQnj173vMQYt++fVv9Wm1tLYWFhZw9e5Zz585x5MgR9u3bR5cuXRg+fDgjRoygf//+uFy6sYaIt/jd737HzJkzef7559m+ffstX6uvr2fTpk2kpaXR0NDA1atXiY6Obvl6aWkpQUFBra6U//rXv2bFihUtj2+stgEsXbqUTZs2tTw+e/Ysf/mXf9nyODY2lv/xP/5HywUppaWlPPHEE7z6atOBvAMHDlBSUsI///M/3/I9v/e9793vH8EDeZCVLxHxENevX2f//v3s37+fkpISgoKCSEhI4PHHHyc2NrZdV6RulLkbJ6M2NDRw8uRJDh8+zKFDh9i3bx+hoaEkJiYybtw4unbt2m7fW0TaX3V1NdXV1bc8V1lZSWNjI6WlpS3P/fu//zvHjx+ntLSULl263PGk9nPnzjFixAiuX7/e8tzLL7/MU0899Y1Vqxu++93vtrry9dJLL7W68lVcXExcXBwAdXV1FBcXs2TJEgIDA/nkk09Ys2YNu3fvZtKkSeTn53PkyBFOnz7N/Pnz2/Tn0h5UvkS8UEVFBbt27WLv3r00NDQwYMAA0tLSGD58eKcttfv5+TF06FCGDh1KXV0dx48f5+DBg+zcuZOMjAwSEhIYP348PXr06JQ8InJ/fvazn/GjH/3ojl+7+Xyrm7399tv8zd/8TcvjmTNn4nK5sNZSW1t7yxWhdXV1rF69mu9+97sAvPjii7z77rstX//Xf/1X/vCHP9zy/iEhIURGRvLee+/dsip2+fJl/vqv/xqAyMhILl26RG1tLQsWLCAjI4O0tDR+8pOfEBcXx5dffsmVK1d47733CAwM5K233uL48eP827/9233+CT04lS8RL1JVVUVGRgbZ2dnU19czatQoJk2a5HjBCQgIID4+nvj4eIqKisjMzCQvL4/c3FyGDx9OamoqvXr1cjSjiNzqzTffvGXlCZoOBa5evZqVK1e26T02bNhAWloaBQUF31j5euqpp+668vXXf/3XLd9/06ZNXLlyhQkTJjBhwgRKSkoYN24c06dPB755zldRUREvv/wyEyZMaDncOW3aNP7qr/6K+vp6ampqKC0tJTo6mtOnTzNw4MA2fZ72ovIl4gXq6+vJzMxk9+7d1NTUMGLECNLS0hwvXXcSFRXFggULmDJlCtnZ2eTk5HDs2DESEhKYNm0aoaGhTkcUaTfr1q3j0qVLjmbo1asXs2fPbpf3unTpEt26dWPAgAFs2LCBoUOHtjq2sbHxlnM8q6ureeqpp1oeZ2dn3/L4bkJCQvjBD37A4sWLuX79Om+//TZ79+6949jCwkJGjx7Nd77zHf7+7/+elStX8q1vfYvvf//7bNmyhbfeeovw8HDWr1/Piy++yMGDB3nyySfb+CfQPlS+RDzcuXPnWLlyJVeuXGHYsGGkpaXxyCOPOB3rnkJDQ5k2bRoTJ05kx44dZGdnc/ToUSZNmsT48ePveDWUiDjr4MGDjBkzho8++oja2tq7jq2trW05zSE2Npa8vLxvjOnTp889v2efPn1azj3r2rUr/v7+BAcHM378ePr168e+fftuGR8TE8PKlStJTk7ms88+IygoiMcff5zw8HCeeeYZsrKySEtL4+OPP2bRokUcOnTolvPJOoP+dhPxUDU1NWzZsoU9e/YQFhbGCy+8wJAhQ5yOdd+Cg4OZOXMmSUlJbNy4kS1btpCbm8usWbMYPny40/FEHkp7rTi5gxt/53z/+99v0/jy8nJCQ0Nv2avrTr71rW/xm9/8ptWvnz9/vuWfX375ZRISEu56VeK5c+eYM2cO0HTRkZ+fH/Hx8S1f//nPf056ejrf//73efPNN5k4cWKr57B1FF3zLeKB8vPz+dWvfsWePXsYO3Ys/+k//SePLF43i4yM5Nlnn+Wll14iODiYTz/9lM8++4zKykqno4kI8Pvf/54+ffqQmJh4z7GXL1+mqqqKnj17UlpayhdffEFYWBjZ2dmUlpZSWlrK+++/jzHmrmXuww8/JDw8vOV/H3zwAT/84Q9veW7z5s23vKZv374UFRXxP//n/2TChAlUVFSQk5ODtZavv/6av/iLv6B79+78t//23/hf/+t/tez91ZlUvkQ8SENDAxs2bOCDDz4gICCAV155hTlz5njVfeQeffRRXnvtNdLS0jh69Ci/+tWvOH78uNOxRHxaYWEhP/zhD/mHf/iHVse88847/OpXTbdszszMpGvXri2nQKSmpvLLX/6S2bNns2LFCn7729/y+uuvs2LFiltWpW73wgsvtJS1K1euEBUVxfTp09mxY0fL89OmffPGOVVVVfzd3/0dvXr14vTp0/zHf/wH6enpt2yv06VLF4wxLF26lLq6ugf9o3kgKl8iHuLatWssXbqUzMxMkpOTef311++66akn8/PzY/Lkybz66quEhITw0Ucf8ac//ekbew6JSMe7fPky8+bNY/LkyTzzzDOtjjtw4ACrV68Gmna+nzZt2i0bN6enp/OLX/yCp59+mldffZXXXnuN0aNHtylDVVUVr7zyCgkJCXzrW99izpw5LFiwgEOHDt1xfJcuXcjMzKRfv34kJibyy1/+krS0tJavL1++nLfeeouMjAzy8/N59tlnuXbtWpuytAed8yXiAQoKCli2bBm1tbUsWrSIkSNHOh2pU/Tu3ZvXXnuN7du3s3v3bk6fPs3ixYu/cRsTEekYGzZs4LXXXmPkyJG8//77AC0Xw2zZsgVjDMaYlnsppqWlcfnyZZYuXcqHH35IXV0dX375JTt37mTFihUcPnyYn/3sZzz66KP8+te/Jjo6mnHjxjFy5EgmT55Menr6NzK8+eab/Pa3vyU1NZVly5bRtWtXFi1axD/+4z+yYsUKamtrOXPmDP369Wt5jbWWkpISzp8/T2RkJG+88QZ/+Zd/SWFhIRcvXuQ3v/kNK1euJCUlhQ0bNjBlyhRGjx7N1q1b23QRwMPSypeIG7PWkpGRwdKlSwkODubVV1/1meJ1g7+/P9OmTeOVV14Bmm5nkpubi7XW4WQi3u2NN95g4cKFLFmyhBUrVrTcBig4OJiXXnqJf/iHf2DkyJGMGDGCBQsWEBISwne/+13Cw8P57ne/y1NPPcV//s//mYULF5KTk8Mbb7zBmTNn+K//9b/yzDPPsHnzZr788ktefPFFiouLWz30N3HiRL744gs++OCDljtjhIaG8vOf/5whQ4bw7LPPcv78eWbOnNnymueff55vfetbDBo0iIMHD/LWW2+Rl5fHlClT2L59O5mZmaSmpgJN54hlZmYyd+5cevfu3cF/qk2MO/8FlpSUZHNycjr0exQUFLTcDkXch+alaY+clStXkpeXx/Dhw3nyyScdPbfLHeaksrKSzz//nJMnT5KQkMDcuXMdvTmuO3CHeRE4duxYy9W5NTU1XnEe5tmzZwFuWVG6X9bae94ztiNUV1ffsps+tM+83DzPd2KMybXWJt3rfXTYUcQN1dXVsWzZMo4fP05qaippaWmO/AXmbkJCQnjhhRfYtm0bO3fu5NKlSyxevLjTLxMX8QUPU7pucOrvrduLl7vRYUcRN1NZWcnSpUs5fvw4c+fOZcqUKSpeN3G5XEydOpXnn3+e0tJS3n33Xc6cOeN0LBGRNlP5EnEjpaWl/O53v+PixYssXryY5ORkpyO5raFDh/Laa68REhLC73//e44cOeJ0JBGRNlH5EnETly9f5re//S0VFRW8+OKL2t29DSIjI3nllVeIjY1l2bJlZGRk6ER8cYz+3fNu7Tm/Kl8ibqCoqIilS5cC8Od//uf079/f4USeIyQkhBdffJH4+Hg2btzI2rVraWxsdDqW+JjAwECqqqqcjiEdqKqqqt0u8NEJ9yIOu3r1asv+Od/5zneIiopyOJHn8ff3Jz09ne7du5OZmUl5eTnp6ek+fyWkdJ6oqCjOnz9PVFQUQUFBBAYG6lxNL2GtpaqqigsXLrTs2P+wVL5EHFRSUsLSpUtpbGxU8XpIxhhmzpxJeHg4a9eu5aOPPuK5554jMDDQ6WjiA8LCwggKCuLKlStcvHgRl0sHltxNfX19ywax9ysgIIBHHnmE7t27t0sWlS8Rh5SVlfH+++9TV1fHSy+9RHR0tNORvMLYsWMJCgriT3/6Ex988AEvvPCCV+y5JO4vODiYvn37au81N+VO86JqLuKA8vJy3n//fWpqanjxxRfp1auX05G8yqhRo0hPT+fcuXP8/ve/1z0hRcStqHyJdLLq6mo++OADKioq+LM/+7NOu52Fr4mPj2fx4sVcvHiRpUuXUllZ6XQkERFA5UukUzU0NPDHP/6RoqIinn32WWJjY52O5NWGDRvGc889x+XLl3n//fepqKhwOpKIiMqXSGex1rJq1SpOnTrFggULGDhwoNORfMKQIUN44YUXKC4u5g9/+IMOQYqI41S+RDrJjh07yMvLIzU1lYSEBKfj+JSBAwfy7LPPcvnyZT788ENqa2udjiQiPkzlS6QTHDhwgG3btjFq1CjS0tKcjuOTBg8eTHp6OufPn+eTTz6hvr7e6Ugi4qNUvkQ62OnTp/niiy949NFHWbBggTZedFBcXBwLFy7k1KlTLFu2jIaGBqcjiYgPUvkS6UClpaX88Y9/pEePHixevBg/Pz+nI/m8hIQE5syZw1dffcWf/vQn3Y9PRDqdNlkV6SB1dXV8+umnNDY28uyzzxIcHOx0JGk2duxYampq2LJlC126dGHOnDlORxIRH6LyJdIBrLWsWbOGixcv8txzz9GjRw+nI8ltnnjiCSoqKsjOziY8PJzx48c7HUlEfITKl0gHyMnJabmy8bHHHnM6jrRi1qxZlJeXs2HDBsLCwoiLi3M6koj4AJUvkXZ27tw51q1bx5AhQ3Rlo5szxrBo0SKuX7/O559/TmhoKP369XM6loh4OZ1wL9KOrl+/zqeffkpYWBiLFi3SlY0eICAggOeee46wsDA+/vhjrl696nQkEfFy9yxfpsm7xph9xpgftzImyBiz3Biz3xjzevNzo4wxGcaYXcaY95vf5xvjRLxFY2Mjy5Yto7q6mmeffZYuXbo4HUnaKCQkhG9/+9sYY1ruuyki0lHasvKVAgQDScASY0zMHcY8DRwG0oCfGmOCAD/gaWvtJKA38Egr40S8wq5duzhz5gzz5s3jkUcecTqO3KfIyEheeOEFrl27pk1YRaRDtaV8jQG2AiOBXUDiXcakAFnAYGvtPmCyMeYEsM9ae+lO4x76E4i4gXPnzrFt2zZGjBjBqFGjnI4jDyg2NpannnqKc+fOsXr1au0BJiIdoi3lKwwoAqKBbc2PWxtjgWM3xlhrPwHigWHGmP6tjRPxZNXV1Xz++eeEhYUxb948nefl4eLj40lNTSUvL4/s7Gyn44iIF2rL1Y5lQBBwlKaVq9OtjIkCDgDpQJkxZjRw0Fpba4zZTtNq1zfG3f5GxpglwBJo+q/QgoKC+/xI90cn17onT5kXay07duygrKyM2bNnc+nSJacjdRhPmZP2MGDAAAoKCtiwYQMNDQ3ExsY6HalVvjQvnkJz4p7caV7aUr5ygXRr7R+NManA+62MiaNpZWwkkA/8C/AJTYcZE4HdNK203T7uFtbad4B3AJKSkuyAAQPu5/M8kM74HnL/PGFe8vLyKCgoYMqUKYwdO9bpOB3OE+akvcTGxvK73/2OXbt28eqrr7r1Rrm+NC+eQnPintxlXtpy2DELCDfGZAN7rbWFdxizHJjWPPY9a20N8BPg74wxu4GL1to9rYwT8UhXr15lzZo19O/fn0mTJjkdR9pZYGAgzz33HC6Xi48++ojq6mqnI4mIl7jnypdtOuP01RuPjTGTgc3AIGvtmeYx1cCi2153Hphx23PfGCfiiRobG/n888/x8/Nj0aJFuFzaMs8bhYeHs3jxYpYuXcrnn3/O888/r3P6ROShPchvjBwgAbjTCpiIT9i9ezeFhYXMnz+fsDBdN+LN+vfvz+zZszlx4gQ7d+50Oo6IeIH7Ll/W2gpr7WFrbV1HBBJxd5cvX2bbtm3ExcURHx/vdBzpBElJSTz++ONs3bqVkydPOh1HRDycjpWI3IeGhgZWrFhBcHAwc+fOdTqOdBJjDPPnzyc6OprPPvuMsrJvXKgtItJmKl8i92H37t1cvHiRefPm0bVrV6fjSCcKDAxk8eLFNDY28umnn2oHfBF5YCpfIm309ddfs337duLj44mLi3M6jjigR48ePPnkkxQWFrJ+/Xqn44iIh1L5EmmDhoYG/vSnPxEcHMycOXOcjiMOGj58OBMmTCAnJ4cDBw44HUdEPJDKl0gb6HCj3GzatGn079+f1atXc+XKFafjiIiHUfkSuYeioiJ27Nihw43SwuVykZ6eTkBAAMuWLaOuThd/i0jbqXyJ3IW1llWrVhEQEMDs2bOdjiNupFu3bjz11FNcvnxZ53+JyH1R+RK5iwMHDnDmzBmmT59OaGio03HEzQwZMoTx48eTm5vLkSNHnI4jIh5C5UukFZWVlWzYsIG+ffsyevRop+OIm5o2bRqxsbGsXLmSkpISp+OIiAdQ+RJpxcaNG6mpqWH+/Pm6n5+0ys/Pj/T0dAA+++wzGhoaHE4kIu5O5UvkDgoKCsjLy2P8+PFER0c7HUfcXEREBAsXLuTChQts3rzZ6Tgi4uZUvkRuU19fz+rVqwkPD2fy5MlOxxEPERcXx5gxY8jMzOTUqVNOxxERN6byJXKb3bt3U1RUxNy5cwkICHA6jniQWbNmERUVxYoVK6isrHQ6joi4KZUvkZuUlJSwc+dO4uPjGTJkiNNxxMMEBATw9NNPU1FRwapVq7DWOh1JRNyQypfITdavX4/L5WLmzJlORxEP1bt3b6ZMmcKxY8fIy8tzOo6IuCGVL5FmJ0+e5KuvviI1NZXu3bs7HUc82IQJE+jfvz/r1q2juLjY6Tgi4mZUvkRounH2unXriIiIICUlxek44uFcLheLFi3CGMPy5ctpbGx0OpKIuBGVLxFgz549FBUVMWvWLPz9/Z2OI14gLCyM+fPnc/78eXbs2OF0HBFxIypf4vMqKirYvn07gwcPZujQoU7HES8yYsQIHn/8cXbs2MGFCxecjiMibkLlS3ze5s2bqaurY9asWdrJXtrd3Llz6datGytWrKCurs7pOCLiBlS+xKcVFhayf/9+xo0bR1RUlNNxxAsFBwezcOFCioqK2LJli9NxRMQNqHyJz7LWsnbtWrp27aqd7KVDDRo0iKSkJLKysjhz5ozTcUTEYSpf4rMOHz7M+fPnmTZtGkFBQU7HES83Y8YMIiIiWLFiBbW1tU7HEREHqXyJT6qvr2fz5s306tWLhIQEp+OIDwgMDOTJJ5+ktLSUDRs2OB1HRByk8iU+KTs7m7KyMmbMmKGT7KXT9O/fn/Hjx5Obm0t+fr7TcUTEISpf4nMqKyvZuXMnQ4YMYeDAgU7HER8zdepUoqKi+OKLL6iurnY6jog4QOVLfM6OHTuora1l+vTpTkcRH+Tv78+iRYu4fv26Dj+K+CiVL/EpxcXF7N27l8TERKKjo52OIz4qJiaGCRMmsH//fk6ePOl0HBHpZCpf4lM2b96Mn58faWlpTkcRH5eWlkZUVBQrV66kpqbG6Tgi0olUvsRnnDt3jqNHjzJx4kS6devmdBzxcf7+/ixcuJCysjI2bdrkdBwR6UQqX+ITrLVs2LCB0NBQxo8f73QcEQD69u1LSkoKOTk5FBQUOB1HRDqJypf4hGPHjnH+/HmmTJlCYGCg03FEWkydOpXIyEi++OILbb4q4iNUvsTrNTY2snXrVqKiorShqridgIAAFi5cSElJie79KOIj7lm+TJN3jTH7jDE/bmVMkDFmuTFmvzHm9ebnHjXGbDbG5BhjftT83BvGmEPGmF3GmLXt+1FE7uzAgQMUFRUxdepUXC7994a4n/79+5OcnEx2djbnzp1zOo6IdLC2/CZKAYKBJGCJMSbmDmOeBg4DacBPjTFBwA+B/26tTQLmGWO6A6HAX1lrJ1lr57THBxC5m/r6erZv305MTAzDhg1zOo5Iq6ZPn05YWBhffPEF9fX1TscRkQ7UlvI1BtgKjAR2AYl3GZMCZAGDgU+Afc1fLwVqaCpff2eM2XNjNUykI+Xk5FBWVsa0adN0GyFxa4GBgcybN4+ioiJ27drldBwR6UD+bRgTBpwFooFtzY/vNKYI6AUcA8KstRubD1m+BXxora0xxhwDNlhrdxpjthpj+ltrz9z8RsaYJcASgNjY2A6/Aujq1asd+v7yYNpjXurq6ti+fTu9evXC5XLparKHpJ+VjhcQEMDAgQPZuXMn4eHhhIeH3/M1mhf3ozlxT+40L20pX2VAEHCUphWu062MiQIOAOnNjwH+N7DDWvspgLX2w5tesxd4FLilfFlr3wHeAUhKSrIDBgxo40d5cJ3xPeT+Pey87Nixg+rqaubNm0efPn3aJ5SP089Kx+vZsyf//u//Tm5uLq+88kqbVmw1L+5Hc+Ke3GVe2nLYMRcYZ63NB1KB/a2MiQOu0HR4Mt8Y8yxw7kbxAjDG/Jsx5nFjjIumc8hOPOwHELmTyspKMjIyGDZsmIqXeJSuXbsye/Zszp8/z969e52OIyIdoC3lKwsIN8ZkA3uttYV3GLMcmNY89j1rbQ3wBrC4+crGXcaYx2laCfsVTWVtvbX2Qrt8CpHb7N69m5qaGqZMmeJ0FJH79vjjjzNo0CA2b95MWVnZvV8gIh7lnocdrbUWePXGY2PMZGAzMOjG+VrW2mpg0W2vS23lLSc9cFqRNigvL2fPnj2MHDlSN88Wj2SMYf78+fzHf/wHq1ev5vnnn9cFIyJe5EE2PcoBEoA7rYCJOG7nzp00Njbq5tni0cLDw5k6dSonTpzgyJEjTscRkXZ03+XLWlthrT1sra3riEAiD6OsrIx9+/aRmJhIRESE03FEHsrYsWOJiYlh3bp1VFVVOR1HRNqJtvsWr7JexzAzAAAgAElEQVRz504AnnjiCYeTiDw8l8vF/PnzqaysZPPmzU7HEZF2ovIlXqO0tJT9+/czevRowsLutB2diOfp3bs3KSkp5ObmcvbsWafjiEg7UPkSr7Fz506MMVr1Eq+TlpZGWFgYq1atoqGhwek4IvKQVL7EK5SUlJCXl8fo0aPp3r2703FE2lVgYCBz587lypUrZGRkOB1HRB6Sypd4hR07dmCMYdIk7WQi3mno0KHExcWxY8cOiouLnY4jIg9B5Us8XnFxMQcOHGDMmDFa9RKvNnv2bFwuF6tXr6ZpC0YR8UQqX+Lxdu7ciZ+fn1a9xOt169aNadOmcerUKQ4dOuR0HBF5QCpf4tGuXr3asurVrVs3p+OIdLikpCRiY2PZsGGD9v4S8VAqX+LRbqx6TZw40ekoIp3C5XIxb948Kisr2bJli9NxROQBqHyJxyopKeHgwYNa9RKf07t3b8aOHUtOTg5XrlxxOo6I3CeVL/FYu3btwuVyMWHCBKejiHS6KVOm0K1bN7KysmhsbHQ6jojcB5Uv8Ujl5eXk5eWRkJCgKxzFJwUFBTFr1iyKi4vZu3ev03FE5D6ofIlH2r17N4CucBSfFhcXR0xMDFu2bOHatWtOxxGRNlL5Eo9z/fp19u3bx8iRIwkPD3c6johjjDGMGzeOhoYG1q9f73QcEWkjlS/xOJmZmTQ0NGjVSwTo3r07TzzxBEeOHCE/P9/pOCLSBipf4lEqKyvZu3cvI0aMoEePHk7HEXELEydOpEePHqxdu5b6+nqn44jIPah8iUfJysqirq6OJ554wukoIm7D39+fuXPnUlxc3HI+pIi4L5Uv8RjV1dXs2bOH4cOH07NnT6fjiLiVgQMHEh8fz65duygpKXE6jojchcqXeIw9e/ZQU1NDamqq01FE3NLMmTNxuVysXbtWN94WcWMqX+IRamtrycrKYsiQIfTq1cvpOCJuqXv37kyePJkTJ07w1VdfOR1HRFqh8iUeYf/+/VRVVekKR5F7GDduHNHR0axbt466ujqn44jIHah8idtraGggIyODfv360a9fP6fjiLg1Pz8/5s6dS1lZGTt27HA6jojcgcqXuL1Dhw5RXl6uVS+RNurfvz+jRo0iIyODoqIip+OIyG1UvsStWWvZvXs3jzzyCIMHD3Y6jojHmD59OgEBAaxZs0Yn34u4GZUvcWtfffUVRUVFTJw4EWOM03FEPEZoaChTp07l9OnTHD161Ok4InITlS9xW9Zadu3aRUREBPHx8U7HEfE4SUlJ9OrVi/Xr11NbW+t0HBFppvIlbqugoIALFy4wYcIEXC79qypyv1wuF3PnzuXatWs6+V7Ejeg3mritXbt2ERoaSkJCgtNRRDxW3759SUhIIDMzUyffi7gJlS9xS1evXuXUqVOkpKTg7+/vdBwRj3bj5HvtfC/iHlS+xC0dOnSIoKAgkpKSnI4i4vG6du3K1KlTOXXqlE6+F3EDKl/idoqLizlz5gxJSUkEBQU5HUfEK9w4+X7Dhg06+V7EYSpf4nYyMzNxuVyMGzfO6SgiXuPGyffl5eU6+V7EYSpf4lYqKirIy8tj0KBBdOvWzek4Il6lb9++jBo1SiffizjsnuXLNHnXGLPPGPPjVsYEGWOWG2P2G2Neb37uUWPMZmNMjjHmR83PRTY/l2eMWdC+H0W8wZ49e6ivr9e+XiIdRCffizivLStfKUAwkAQsMcbE3GHM08BhIA34qTEmCPgh8N+ttUnAPGNMd+B14EPgz4C3Hz6+eJPa2lr27t3LY489RlhYmNNxRLxSaGgoU6ZM4dSpU3z55ZdOxxHxSW0pX2OArcBIYBeQeJcxKUAWMBj4BNjX/PVSoOamcXFAsTGm68OEF++Sl5dHVVUVEyZMcDqKiFdLTk4mOjqa9evXU1dX53QcEZ/Tlg2UwoCzQDSwrfnxncYUAb2AY0CYtXZj8yHLt4APrbU1xpgb40qAwubXVdz8RsaYJcASgNjYWAoKCh7gY7Xd1atXO/T9pW0aGxvZuXMnPXv2pLGxUfPihjQn7ulB5yUxMZH169ezatUqEhPv9N/U8qD0s+Ke3Gle2lK+yoAg4ChNK1enWxkTBRwA0psfA/xvYIe19tObxg2iaQXtb24a18Ja+w7wDkBSUpIdMGBAGz/Kg+uM7yF3d/jwYa5fv868efNa5kPz4n40J+7pQeZlwIABFBYWcuTIESZPnkxkZGT7B/Nh+llxT+4yL2057JgLjLPW5gOpwP5WxsQBV2g6PJlvjHkWOHdT8boxLrb5+3az1lZ8453E51hrycjIoEePHjz22GNOxxHxGTNmzMDPz49169Y5HUXEp7SlfGUB4caYbGCvtbbwDmOWA9Oax75nra0B3gAWG2N2Nf/vcZpWtL5P0+HLf2qPDyCer6CggIsXLzJ+/HiMMU7HEfEZ3bp1Y/LkyZw4cYLjx487HUfEZ9zzsKNtuhb51RuPjTGTgc3AIGvtmeYx1cCi216X2spbTnvgtOKVMjIy6Nq1K6NGjXI6iojPGTduHPv372fdunUMHDhQ91IV6QQPsslqDpBA0wnzIg/l8uXL5OfnM3bsWP2lL+IAPz8/5syZQ0lJCRkZGU7HEfEJ912+rLUV1trD1lpdnywPLTMzE39/f91AW8RBAwcOJC4ujp07d1JaWup0HBGvp9sLiWOuX7/OoUOHSEhIICQkxOk4Ij5t5syZGGPYsGGD01FEvJ7Klzhmz549NDQ0MH78eKejiPi8sLAwnnjiCY4dO8bJkyedjiPi1VS+xBG1tbXk5OQwbNgw7S8k4ibGjx9PZGQk69ato6Ghwek4Il5L5UscceDAAaqqqrTqJeJG/P39mT17NkVFRWRnZzsdR8RrqXxJp2tsbCQrK4vY2Fj69u3rdBwRucmQIUMYOnQo27dv59q1a07HEfFKKl/S6Y4fP05xcbE2VRVxU7NmzaKhoYGNGzc6HUXEK6l8SafLzMwkLCyM4cOHOx1FRO4gMjKSiRMncujQIc6cOeN0HBGvo/IlnerChQucPXuWlJQUXC796yfiriZNmkRYWBhr166lsbHR6TgiXkW//aRTZWZmEhQURGJiotNRROQuAgICmDVrFl9//TU5OTlOxxHxKipf0mlKS0s5evQoY8aMISgoyOk4InIPw4YNY+DAgWzZsoWKigqn44h4DZUv6TQ3Ll0fO3asw0lEpC2MMcyZM4e6ujo2bdrkdBwRr6HyJZ2ipqaG/fv3Ex8fT1hYmNNxRKSNoqKiSElJIS8vj/PnzzsdR8QrqHxJp9i/fz81NTWkpKQ4HUVE7lNqaiqhoaGsWbNGJ9+LtAOVL+lwjY2NZGdn07dvX2JjY52OIyL3KSgoiJkzZ3Lx4kX279/vdBwRj6fyJR3uq6++orS0VKteIh5sxIgR9OvXj82bN1NZWel0HBGPpvIlHS4rK4vw8HCGDRvmdBQReUDGGObOnUt1dTVbt251Oo6IR1P5kg5VWFjI2bNnGTt2rDZVFfFwjzzyCMnJyeTk5HDx4kWn44h4LP02lA6VlZVFYGAgo0ePdjqKiLSDKVOm0LVrV9asWYO11uk4Ih5J5Us6THl5OUeOHCExMVGbqop4ieDgYKZPn8758+fJy8tzOo6IR1L5kg6zZ88erLWMGzfO6Sgi0o5GjRpFnz592LRpE9XV1U7HEfE4Kl/SIWpra8nNzWXYsGFEREQ4HUdE2tGNk++rqqp08r3IA1D5kg5x4MABqqurtb2EiJfq3bs3Y8aMYe/evVy6dMnpOCIeReVL2p21luzsbGJiYujbt6/TcUSkg0ydOpXg4GDWrl2rk+9F7oPKl7S7/Px8rl69yrhx4zDGOB1HRDpIly5dmD59OmfPnuXgwYNOxxHxGCpf0u6ys7MJDQ0lPj7e6Sgi0sESExOJjY1l48aNOvlepI1UvqRdXb58mZMnT5KcnIyfn5/TcUSkg904+b6iooJt27Y5HUfEI6h8SbvKzs7G39+fpKQkp6OISCeJiYlhzJgx7NmzRyffi7SBype0m8rKSg4ePMjIkSMJCQlxOo6IdKJp06YRHBysne9F2kDlS9pNbm4u9fX12lRVxAfdOPn+3LlzOvle5B5UvqRdNDQ0sHfvXgYOHEh0dLTTcUTEATr5XqRtVL6kXRw9epRr165pU1URH2aMYd68eVRUVGjne5G7UPmSh2atJSsrix49ejB48GCn44iIg3r37k1SUpJ2vhe5i3uWL9PkXWPMPmPMj1sZE2SMWW6M2W+Meb35OT9jzE+NMbtuGveGMeaQMWaXMWZt+30McdL58+cpLCzUpqoiAjTtfN+lSxedfC/SirasfKUAwUASsMQYE3OHMU8Dh4E04KfGmCDAD9gN+N80LhT4K2vtJGvtnIcJLu4jKyuL4OBgRo0a5XQUEXEDXbp0YcaMGZw7d468vDyn44i4nbaUrzHAVmAksAtIvMuYFCALGGytrbXWrr5tXCjwd8aYPcaYHz14bHEXZWVlHDt2jNGjRxMYGOh0HBFxE6NGjaJfv35s3LiRyspKp+OIuJW2lK8woAiIBrY1P25tjAWOtTKG5q/9o7V2LJBqjOl/v4HFvezZsweAsWPHOpxERNzJjZ3vq6ur2bx5s9NxRNyK/72HUAYEAUdpWuE63cqYKOAAkN78+BustR/e9HAv8Chw5uYxxpglwBKA2NhYCgoK2hDxwV29erVD39+b1dXVkZOTQ79+/SgpKaGkpKTd3lvz4n40J+7J3edl+PDh7Nu3j169etGzZ0+n43QKd58TX+VO89KW8pULpFtr/2iMSQXeb2VMHE0rYyOB/Du9kTHm34D/Axyh6Ryyf7l9jLX2HeAdgKSkJDtgwIA2RHw4nfE9vFFOTg61tbVMnTqVfv36tfv7a17cj+bEPbnzvPTu3Ztz586xb98+XnvtNVwu37jI3p3nxJe5y7y05acgCwg3xmQDe621hXcYsxyY1jz2PWttTSvv9b+BX9FU1tZbay88QGZxA9ZasrOziYmJoW/fvk7HERE3FRQUxKxZs7h06RI5OTlOxxFxC/dc+bJN1wm/euOxMWYysBkYZK090zymGljUyutTbvrnr4BJD5lZ3MDJkycpKipi0aJF2l5CRO4qLi6OgQMHsmXLFuLi4ggNDXU6koijHmT9NwdIAO60AiY+Ijs7m9DQUOLj452OIiJu7sbJ9/X19WzYsMHpOCKOu+/yZa2tsNYettbWdUQgcX9XrlwhPz+f5ORk/Pz8nI4jIh6gR48eTJw4kUOHDnH69J2u2xLxHb5x5qO0q+zsbPz8/BgzZozTUUTEgzzxxBNERESwevVq6uvrnY4j4hiVL7kvVVVVHDx4kMcff5yuXbs6HUdEPIi/vz/z5s3j6tWr7N692+k4Io5R+ZL7sm/fPurq6khJSbn3YBGR2wwaNIj4+Hh27txJcXGx03FEHKHyJW3W0NDAnj17ePTRR3nkkUecjiMiHmrWrFn4+/vrxtvis1S+pM2+/PJLysvLGTdunNNRRMSDdevWjalTp3Ly5EmOHDnidByRTqfyJW2WlZVFREQEQ4cOdTqKiHi4pKQkevfuzfr166murnY6jkinUvmSNrlw4QLnz59n3Lhx2lRVRB6ay+Vi/vz5VFRUsGXLFqfjiHQqlS9pk+zsbIKCgkhISHA6ioh4iZiYGJKTk9m7dy/nz593Oo5Ip1H5knsqLy/nyJEjJCYmEhQU5HQcEfEiU6dOpVu3bqxatYqGhgan44h0CpUvuae9e/dirWXs2LFORxERLxMUFMScOXP4+uuvycrKcjqOSKdQ+ZK7qqurIzc3l8cee4yIiAin44iIFxo+fDiPPfYY27Zto6SkxOk4Ih1O5Uvu6uDBg1RVVWlTVRHpUHPmzMHlcmnvL/EJKl/SKmst2dnZ9OrVi379+jkdR0S8WFhYGFOnTiU/P197f4nXU/mSVp06dYorV66QkpKi7SVEpMMlJycTExPDunXrqKqqcjqOSIdR+ZJWZWdn07VrV+Lj452OIiI+4MbeX5WVlWzatMnpOCIdRuVL7qioqIgTJ06QnJyMv7+/03FExEf07t2blJQU9u3bR0FBgdNxRDqEypfcUVZWFn5+fiQlJTkdRUR8zJQpU4iIiGDlypXU19c7HUek3al8yTdUVlZy4MABRo4cSdeuXZ2OIyI+JiAggPnz51NcXMz27dudjiPS7lS+5Btyc3Opr6/X9hIi4piBAweSkJBARkYGly5dcjqOSLtS+ZJbNDQ0sHfvXgYOHEh0dLTTcUTEh82cOZMuXbqwcuVKGhsbnY4j0m5UvuQWR44c4dq1a1r1EhHHdenShTlz5lBYWEh2drbTcUTajcqXtLDWkpWVRVRUFIMHD3Y6jogIcXFxDB06lC1btujWQ+I1VL6kxdmzZ7l48SLjxo3Tpqoi4haMMcybNw+Xy8WqVat06yHxCipf0iIrK4suXbowatQop6OIiLTo3r07M2bM4NSpU+Tl5TkdR+ShqXwJACUlJXz55ZeMGTOGgIAAp+OIiNxizJgxDBgwgPXr11NeXu50HJGHovIlQNOthFwuF8nJyU5HERH5BmMMCxYsoLGxUYcfxeOpfAnV1dXs37+f+Ph4unfv7nQcEZE7ioyMZOrUqZw4cYJDhw45HUfkgal8Cbm5udTW1jJ+/Hino4iI3NXYsWPp27cva9eu5fr1607HEXkgKl8+rqGhgT179jBgwAB69+7tdBwRkbtyuVwsXLiQuro6Vq9ercOP4pFUvnzckSNHKC8v16qXiHiMqKgopkyZwpdffsnRo0edjiNy31S+fNjNm6oOGTLE6TgiIm02fvx4YmJiWLNmDRUVFU7HEbkvKl8+7MyZM1y8eJGUlBRtqioiHsXlcvHkk09SU1Ojw4/icVS+fFhmZiYhISGMHDnS6SgiIvctOjqaKVOmcOzYMQ4fPux0HJE2u2f5Mk3eNcbsM8b8uJUxQcaY5caY/caY15uf8zPG/NQYs+umcZHGmM3GmDxjzIL2+xhyv4qKijh+/DjJycnaVFVEPNb48ePp06cPa9as4dq1a07HEWmTtqx8pQDBQBKwxBgTc4cxTwOHgTTgp8aYIMAP2A343zTudeBD4M+Atx88tjysrKws/Pz8tKmqiHg0l8vFU089RX19PStXrtThR/EIbSlfY4CtwEhgF5B4lzEpQBYw2Fpba61d3cq4OKDYGNP1QYPLg6uoqODAgQOMHDmSrl01BSLi2Xr06MH06dM5ceIE+/fvdzqOyD21pXyFAUVANLCt+XFrYyxwrJUxN48rAQrvMk46UE5ODvX19dpeQkS8xtixY1vu/VhaWup0HJG78r/3EMqAIOAoTStXp1sZEwUcANKbH7f2XoNoWkH7mzuNM8YsAZYAxMbGUlBQ0IaID+7q1asd+v7upr6+nqysLGJjY6moqHDbS7R9bV48gebEPWle/r8xY8Zw/vx5PvnkE2bOnOnYVdyaE/fkTvPSlvKVC6Rba/9ojEkF3m9lTBxNK2Mjgfy7vFcscBzoZq39xm9+a+07wDsASUlJdsCAAW2I+HA643u4i5ycHKqrq5k+fbrbf253z+eLNCfuSfPy/9XU1LBq1Sq+/vprUlJSHMuhOXFP7jIvbTnsmAWEG2Oygb3W2sI7jFkOTGse+561tqaV93oH+D5NJe2f7j+uPIzGxkYyMzOJiYmhf//+TscREWl3o0ePZujQoWzatInLly87HUfkju5ZvmyTV62146y1f2+MmWyMqTfG9L9pTLW1dpG1dqy19v/c9vqUm/75qrV2mrU22Vq7sn0/itzLV199RXFxMRMmTNCmqiLilYwxLFiwgKCgIJYvX059fb3TkUS+4UE2Wc0BEmg6YV48hLWW3bt3ExERwfDhw52OIyLSYUJDQ1m4cCGXLl1i69atTscR+Yb7Ll/W2gpr7WFrbV1HBJKOcfbsWS5cuMD48eNxuXRjAxHxbo899hijR48mIyOjwy/cErlf+i3sIzIyMggJCSEhIcHpKCIinWLWrFlERkayYsUKqqurnY4j0kLlywdcuXJFtxISEZ8TGBjIokWLKC8vZ+3atU7HEWmh8uUDMjIy8Pf3Z+zYsU5HERHpVH369CE1NZWDBw9y6NAhp+OIACpfXq+8vJyDBw+SmJhISEiI03FERDpdamoqffr0YfXq1ZSUlDgdR0Tly9tlZ2djrdWthETEZ7lcLtLT0wH47LPPaGhocDiR+DqVLy9WXV1Nbm4ucXFxREREOB1HRMQx4eHhLFiwgAsXLmj7CXGcypcX27t3LzU1NUycONHpKCIijouPjycxMZHdu3dz6tQpp+OID1P58lJ1dXVkZWUxePBgevfu7XQcERG3MHv2bKKioli+fDkVFd+4vbBIp1D58lL79++nsrKSSZMmOR1FRMRtBAYGkp6eTlVVFStWrMBa63Qk8UEqX16ooaGBjIwM+vbtqxtoi4jcplevXsycOZP8/HyysrKcjiM+SOXLCx06dIiysjKeeOIJp6OIiLil5ORkhg0bxqZNmzh37pzTccTHqHx5mRs30H7kkUcYPHiw03FERNySMYYnn3yS7t27s2zZMiorK52OJD5E5cvLfPnllxQVFTFp0iSMMU7HERFxW8HBwTzzzDNUVFSwfPlynf8lnUbly4tYa9m5cyeRkZHExcU5HUdExO3FxMQwa9Ys8vPz2bVrl9NxxEeofHmRU6dOcfHiRSZOnIjLpakVEWmLpKQkRowYwdatWzl9+rTTccQH6De0F9m1axfdunVj5MiRTkcREfEYxhjmz59PZGQkn332GdevX3c6kng5lS8vcfbsWQoKChg/fjz+/v5OxxER8ShBQUEsXryYmpoali1bpvs/SodS+fIS27dvp2vXriQlJTkdRUTEI0VHR7NgwQLOnDnDxo0bnY4jXkzlywucO3eOU6dOMWHCBAICApyOIyLisUaOHMm4cePIzs7m4MGDTscRL6Xy5QW2b99OSEiIVr1ERNrBjBkz6N+/PytXruTixYtOxxEvpPLl4c6fP8/JkycZP348gYGBTscREfF4fn5+PPPMM4SEhPDJJ59oA1ZpdypfHm779u106dKFsWPHOh1FRMRrdO3alcWLF3P9+nWWLVtGY2Oj05HEi6h8ebALFy6Qn5+vVS8RkQ4QGxvLvHnzOH36NJs2bXI6jngR7UngwbTqJSLSsRITEyksLCQzM5Po6GgSEhKcjiReQCtfHqqwsJATJ06QkpJCUFCQ03FERLzW7NmzefTRR1m5ciVnz551Oo54AZUvD7V9+3aCg4MZN26c01FERLzajRPwIyIi+OSTTygpKXE6kng4lS8PVFhYyPHjx7XqJSLSSbp06cLzzz9PY2MjH330ETU1NU5HEg+m8uWBtmzZQpcuXbTqJSLSiXr06MEzzzxDUVERn332ma6AlAem8uVhCgoKOHnyJJMmTSI4ONjpOCIiPmXgwIHMnTuXEydOsGHDBqfjiIfS1Y4exFrLli1b6NatG8nJyU7HERHxSUlJSRQVFZGdnU14eDgpKSlORxIPo5UvD3LixAnOnTtHamqq7uEoIuKgmTNnMmzYMNavX8/Ro0edjiMeRuXLQ9xY9YqIiCAxMdHpOCIiPs3lcvH000/Tt29fPv/8c86cOeN0JPEgKl8e4siRI3z99dekpaXh5+fndBwREZ8XEBDAc889R3h4OB9//DFXrlxxOpJ4iHuWL9PkXWPMPmPMj1sZE2SMWW6M2W+Meb35uUhjzGZjTJ4xZkHzc28bY/YaY3YZY/5v+34U79XQ0MDWrVuJjo7m8ccfdzqOiIg0CwkJ4dvf/jZ+fn588MEHXLt2zelI4gHasvKVAgQDScASY0zMHcY8DRwG0oCfGmOCgNeBD4E/A95uHhcKPGOtnWSt/fOHzO4zDhw4QHFxMVOnTsUY43QcERG5SUREBC+88AKVlZV8+OGH1NbWOh1J3FxbytcYYCswEtgF3OmEoxtjUoAsYPBNz8UBxcaYrjSVr381xuQaY5Y8fHzvV19fz/bt2+nTpw9Dhw51Oo6IiNxBTEwMixcv5vLly2zevFkFTO6qLeUrDCgCooFtzY9bG2OBY82PbzxXAhQ2P84FfkBTSXvDGKONqu4hOzub8vJyrXqJiLi5wYMHk56ezpUrV/j000+pr693OpK4qbbs81UGBAFHaVrNOt3KmCjgAJDe/LgMGETTatnfAGXW2n++8QJjzFdAL6Dg5jdqXhFbAhAbG0tBwS1fbndXr17t0Pd/GNXV1S2rXsaYDv+zcCfuPC++SnPinjQv7iUkJISRI0dy4MABfv/73zN58mRcLl3b5g7c6WelLeUrF0i31v7RGJMKvN/KmDiaVsZGAvnNz8UCx4Fu1toKY8xnNJ0Ldh0YCFy8/Y2ste8A7wAkJSXZAQMG3OdHun+d8T0exOrVq6mvr2fhwoX07NnT6Tidzl3nxZdpTtyT5sX99OrVi/Xr13Pw4EGefPJJHblwE+7ys9KWOp4FhBtjsoG91trCO4xZDkxrHvuetbaGpgL1fZoK2T81j/slsLp53L80j5M7uHz5Mrm5uSQlJflk8RIR8WQpKSmkpaVx4MAB1q5di7XW6UjiRu658mWb/o159cZjY8xkYDMwyFp7pnlMNbDottddpamQ3fzcbkB3g26DjRs3EhgYSFpamtNRRETkAaSmplJTU0NmZiYul4tZs2ZpBUyAB7u3Yw6QQNNJ9NIB8vPzyc/PZ8aMGYSEhDgdR0REHoAxhhkzZtDY2Eh2djbWWmbPnq0CJvdfvqy1FTTt6SUdoLGxkQ0bNhAREcHYsWOdjiMiIg/BGMOsWbMAWgrYnDlzVMB83IOsfEkH2rdvH1euXGHx4sX4+2t6REQ83Y0CZowhKysLay1z585VAfNh+u3uRqqrq9m6dSv9+/dn2LBhTscREZF2Yoxh5syZGGPIzMwEUAHzYSpfbmT79u1UVoTYX9MAAA6lSURBVFa2/ICKiIj3uHEOmDGGjIwM6urqWLBgAX5+fk5Hk06m8uUmLl36f+3da3BV9X7G8e8/YbOTkLvpCdckSCYqBJQkmMNFgkA1RRG8QGUGYkeBY89MX3TGXmY65/RMOz190Rdn2k5nPCkiqKhVKVKCOkUuQsBASAUDBw0gxIAgMTcM5LID/77YkQEMZEey11qwns/MnmRt/mE/4Ud2nll7rbXPsnfvXgoKChg5sq+3zxQRkdudMYa5c+cydOhQduzYQUdHB8888wyBQMDtaOIgXXbXA6y1VFRUEB8fz9y5c92OIyIiUWSMoaSkhHnz5lFXV8cbb7xBZ2en27HEQSpfHlBTU8Pp06d55JFHiI+PdzuOiIg4YMqUKTz99NOcOnWKNWvW0N7e7nYkcYjKl8va29vZunUrOTk5TJo0ye04IiLioPz8fJYsWUJzczOrV6+mubnZ7UjiAJUvl23ZsoXu7m4ee+wxHWQvIuJDubm5lJWV0dnZyapVq6ivr3c7kkSZypeLTpw4weeff86MGTPIyMhwO46IiLhk9OjRLF++nISEBF5//XUOHjzodiSJIpUvl/T09LB582bS0tKYMWOG23FERMRl6enpvPDCC4wZM4b333+fbdu26Q2571AqXy6prKykqamJefPm6RRjEREBID4+nqVLlzJ58mR27drF+vXrCYVCbseSQabrfLng9OnT7Ny5k4kTJ5Kbm+t2HBER8ZDY2Fjmz5/PXXfdxccff0xTUxOLFy8mLS3N7WgySLTny2GhUIgNGzaQlJTEvHnz3I4jIiIeZIxh+vTpLFmyhNbWVsrLy6mrq3M7lgwSlS+HbdmyhaamJhYsWEBcXJzbcURExMPy8vJYuXIlqampvPXWW2zbto3Lly+7HUtukcqXg44fP051dTXFxcXcfffdbscREZHbQFpaGs8//zwPPPAAu3btYt26dVy4cMHtWHILVL4c0tHRwcaNG8nIyGDOnDluxxERkdtIIBBgwYIFzJ8/n/r6el5++WWOHTvmdiz5iVS+HPLBBx9w4cIFnnrqKZ3dKCIiP0lBQQHLly8nPj6edevW8eGHH+psyNuQypcDamtrOXToECUlJYwYMcLtOCIichsbPnw4K1asoLi4mH379lFeXs6ZM2fcjiUDoPIVZefOnWPTpk2MHj1aF1MVEZFBEQgEKC0tZenSpXR1dbFq1Sp27tzJpUuX3I4mEVD5iqKOjg7efvttgsEgixYtIiZG/9wiIjJ4xo0bx4svvsi9997L9u3bKS8vp6Ghwe1Y0g+1gSi5fPky69evp62tjcWLF5OcnOx2JBERuQMlJCSwaNEinn32WTo7O1m9ejWbN2+ms7PT7WhyA7rCfZRs3bqV48ePM3/+fMaMGeN2HBERucPdc8895OTksH37dvbt28cXX3xBaWkp48ePxxjjdjy5ivZ8RcGhQ4fYs2cPRUVFFBQUuB1HRER8IhgMUlpayvLly0lMTOS9997j1Vdf5dSpU25Hk6uofA2ys2fPsnHjRrKysigtLXU7joiI+NDIkSNZsWIFjz/+OM3NzbzyyiusX7+e1tZWt6MJetlxUDU3N/Pmm29eef09NjbW7UgiIuJTMTExFBYWkp+fz+7du/n00085cuQIxcXFTJs2jWHDhrkd0bdUvgZJS0sLa9eupaenh+eee47ExES3I4mIiBAMBpk9ezZFRUVs27aNPXv2UF1dTWFhIdOmTSMpKcntiL6j8jUIWltbWbt2LaFQiLKyMjIzM92OJCIico3k5GQWLlzIjBkzqKysZO/evVRXV1NQUMD06dNJSUlxO6JvqHzdoh+KV1dXF2VlZQwfPtztSCIiIjeUkZHBwoULmTlzJpWVldTU1FBTU8N9993Hgw8+yJgxY3R2ZJSpfN2CtrY2XnvtNTo6OigrK9NbB4mIyG0jPT2dJ554gpKSEqqqqjhw4ACHDx8mMzOTKVOmMHHiRIYOHep2zDuSytdPdObMGd555x06OjpYtmwZI0eOdDuSiIjIgKWkpPDoo48ye/Zsamtr2bdvHxUVFWzZsoXx48czadIksrOztTdsEKl8DZC1lpqaGj766CMSEhJYtmwZo0aNcjuWiIjILQkEAhQUFDB58mQaGhqoqanh8OHDfPbZZyQnJ5Ofn8/EiRPJzMxUEbtFKl8D0N3dTUVFBbW1tYwbN44nn3xSp+qKiMgdxRhDVlYWWVlZhEIhvvzyS2pra6mqqmLPnj2kpqaSm5tLXl4eOTk5BAIBtyPfdlS+InTu3DneffddmpqaePjhh3nooYfU/EVE5I4WCATIz88nPz+fixcvcuTIEY4ePcrBgwfZv38/Q4YMYezYsYwdO5asrCxGjBhBTIyu396ffsuXCTeMcqAQqLDW/rqPNUHgbSAHeNla+3tjTDrwLnAX8Ctr7SZjzN3AG0A8sNJaWz1o30mUtLS0sHv3bg4cOEBcXBzLli1j7NixbscSERFxVEJCAoWFhRQWFtLT08PJkyepq6vj2LFjHD16FAiXtdGjR5OVlcWoUaPIzMwkKSlJOyuuE8mer58DcUAR8I0x5mVr7TfXrXkKOAT8GfCVMWYN8AvgTWAv8B6wCfgb4FeEy9c/Ap59/53GxkYqKyupra0lJiaG+++/n1mzZulidCIi4ntDhgwhNzeX3NxcAM6fP09DQwP19fU0NDTwySefXFkbHx9PZmYmmZmZZGRkkJ6eTlpaGikpKb7dSxZJ+SoEtgOTgEpgMnB9+SoEPiBc1KqA3N77/ppwaWs2xgwD7gf+AvhTIH0Q8g8Kay2tra00NjbS2NjI119/TV1dHYFAgOLiYqZOnUpycrLbMUVERDwpOTmZCRMmMGHCBAC6uro4e/Ys33777ZVbTU0NPT09V74mJiaGlJQUUlNTSUxMZNiwYSQmJpKYmEhCQgJxcXEEg8ErH4cOHXrH7EGLpHylAF8DPwN29G73teY7YDhwpHf7h/taCJe1FMAQfjPvk8ClW0o+CE6cOEFFRQXnz5+/5j9EcnIyM2fOpLi4mISEBBcTioiI3H6CwSDZ2dlkZ2dfue/y5ct8//33tLS00NzcTEtLC62trbS1tXHq1Cna29sJhUI3/XtjY2MZMmTINTdjDDExMdd8vL6kGWPIyckhJycnGt/ugEVSvtqAIPAHwnuzTtxgTQZwEHi6d7sNGEd4b9lLvdvNQF7vfX2WL2PMSmAlwKhRozh58mTE38xANTU1ERsbS15e3pX2nZqaeuWicufOnYvaY8vNNTU1uR1BrqOZeJPm4j2aSf/S09NJT//xC2ChUIiOjg46OzsJhUKEQiG6u7vp7u6mp6eHS5cu/ehmrf3RDbjy8YfPL168GNVOMRCRlK8a4Glr7bvGmJnA2husGU94z9gk4FjvfaOAOiDJWnvBGHOs9zFHAuf7ejBrbTnhA/wpKiqy0WypOTk5ZGRkeKYJy7U0F+/RTLxJc/EezcR7Tp486Zm5RHKkWxWQaozZC1T3cbA9wAZgTu/aNdbaLsIF6i8JF7J/7l33O+A/gP8GfnNLyUVERERuQ/3u+bLh/XbLf9g2xpQAW4Fx1tr63jWdwJPXfV0T4UJ29X1fAVNvPbaIiIjI7emnnOO5H3iAH5/xKCIiIiL9GPAV7q21Fwhf00tEREREBsifVzcTERERcYnKl4iIiIiDVL5EREREHKTyJSIiIuIglS8RERERB6l8iYiIiDhI5UtERETEQSpfIiIiIg5S+RIRERFxkAm/daM3GWMagfooP0wG8F2UH0MGTnPxHs3EmzQX79FMvMmJuWRba/+ov0WeLl9OMMbst9YWuZ1DrqW5eI9m4k2ai/doJt7kpbnoZUcRERERB6l8iYiIiDhI5QvK3Q4gfdJcvEcz8SbNxXs0E2/yzFx8f8yXiIiIiJO050tERETEQb4pXybsP40x/2eM+YcbrAkaYzYYYz4zxvzC6Yx+E+FM0o0xm4wxnxpjXjHGGKdz+kkkM7lq7UvGmCqnsvlZpHMxxvybMWafMcYzL6/cqQb4/LXLGPOeMSbe6Zx+Y4xZYoz5zhgTd4M/j/g5Lpp8U76AnwNxQBGw0hgzso81TwGHgFnAb40xQefi+VIkM/kl8Lq1diqQCkx0MJ8fRTITjDETgHucDOZz/c7FGDMdaLPWPhjeNBMczug3kfysLAE2W2sfAtoBzST6vge+vMmfR/QcF21+Kl+FwHZgElAJTL7Jmp8DVUCuY+n8KZKZbAX+t/fzbqDFmWi+1e9MjDEB4J+Al5yN5muR/Kz8MRAwxmwDvrHWHnYwnx9FMpOjwL8bY+qARGvtfgfz+ZK1tgII3WRJJHOLOj+VrxTCV7b9GbCjd/tGayxw5AZrZPD0OxNr7afW2lZjzC+BI9baBmcj+k4kPye/Bv7VWtvmYC6/i2Quw4E4a+1s4D5jzBTn4vlSJDM5Csyy1uYBjcaYuc7FkxuIZG5R56fy1QYEga+ApN7tvtZkAAdvskYGTyQzwRjzt0DQWuva6/M+EslM/gT4e2PMDmC8MeavnIvnW5HMpZ3wLxOAT4DxjiTzr0hm8lvghz2Q/wNMdyaa3EREv3eizU/lqwYottYeA2YCn91gzXigkfAuyWPOxfOlfmdijJkKZFlrf+d0OJ/qdybW2iJr7Sxr7SzgD9baf3E4ox9F8vxVBRT3fl7AzY97kVsXyUwApvV+nAocdyKY3FSkc4sqP5WvKiDVGLMXqLbWftPHmg3AnN61a6y1XU4G9KFIZvLnQIkxprL39oizEX0nkpmI8yKZy/vAaGPMHqDVWqszUaMrkpn8Bvi73pmMAf7LwXzSN088x/n2IqvGmBLCB3OPs9bWu51HNBMv0ky8SXPxHs3Ee4wx2YT3Ns6x1n7idp6r+bl8DQPGAl9aa292ZoQ4RDPxHs3EmzQX79FMvKf3zOx7gBPW2gtu57mab8uXiIiIiBv8dMyXiIiIiOtUvkREREQcpPIlIiIi4iCVLxEREREHqXyJiIiIOEjlS0RERMRB/w8ffhoFe+/9SQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ps = np.linspace(0, 1, 100)\n",
    "Ls = [np.prod(stats.bernoulli(prob).pmf(coin_result))\n",
    "      for prob in ps]\n",
    "\n",
    "fig = plt.figure(figsize=(10, 6))\n",
    "ax = fig.add_subplot(111)\n",
    "ax.plot(ps, Ls, label='尤度関数', color='gray')\n",
    "ax.legend(fontsize=16)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.336101Z",
     "start_time": "2018-08-20T18:22:38.330723Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-3.365"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prob = 0.4\n",
    "rv = stats.bernoulli(prob)\n",
    "mll = np.sum(np.log(rv.pmf([0, 1, 0, 0, 1])))\n",
    "mll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.342469Z",
     "start_time": "2018-08-20T18:22:38.337446Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-76.325"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rv = stats.norm(y_hat, np.sqrt(unexp_var / n))\n",
    "mll = np.sum(np.log(rv.pdf(y)))\n",
    "mll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.346776Z",
     "start_time": "2018-08-20T18:22:38.343621Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "156.650"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "aic = -2 * mll + 2 * (p+1)\n",
    "aic"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.351446Z",
     "start_time": "2018-08-20T18:22:38.348026Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "158.642"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "bic = -2 * mll + np.log(n) * (p+1) \n",
    "bic"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## モデルの妥当性"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.369998Z",
     "start_time": "2018-08-20T18:22:38.352673Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.756</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.727</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   26.35</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Tue, 21 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>6.19e-06</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>03:22:38</td>     <th>  Log-Likelihood:    </th> <td> -73.497</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   153.0</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    17</td>      <th>  BIC:               </th> <td>   156.0</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     2</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   -1.8709</td> <td>   11.635</td> <td>   -0.161</td> <td> 0.874</td> <td>  -26.420</td> <td>   22.678</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>      <td>    6.4289</td> <td>    0.956</td> <td>    6.725</td> <td> 0.000</td> <td>    4.412</td> <td>    8.446</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>睡眠時間</th>      <td>    4.1917</td> <td>    1.778</td> <td>    2.357</td> <td> 0.031</td> <td>    0.440</td> <td>    7.943</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 2.073</td> <th>  Durbin-Watson:     </th> <td>   1.508</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.355</td> <th>  Jarque-Bera (JB):  </th> <td>   1.716</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.660</td> <th>  Prob(JB):          </th> <td>   0.424</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.437</td> <th>  Cond. No.          </th> <td>    38.0</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.756\n",
       "Model:                            OLS   Adj. R-squared:                  0.727\n",
       "Method:                 Least Squares   F-statistic:                     26.35\n",
       "Date:                Tue, 21 Aug 2018   Prob (F-statistic):           6.19e-06\n",
       "Time:                        03:22:38   Log-Likelihood:                -73.497\n",
       "No. Observations:                  20   AIC:                             153.0\n",
       "Df Residuals:                      17   BIC:                             156.0\n",
       "Df Model:                           2                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     -1.8709     11.635     -0.161      0.874     -26.420      22.678\n",
       "小テスト           6.4289      0.956      6.725      0.000       4.412       8.446\n",
       "睡眠時間           4.1917      1.778      2.357      0.031       0.440       7.943\n",
       "==============================================================================\n",
       "Omnibus:                        2.073   Durbin-Watson:                   1.508\n",
       "Prob(Omnibus):                  0.355   Jarque-Bera (JB):                1.716\n",
       "Skew:                           0.660   Prob(JB):                        0.424\n",
       "Kurtosis:                       2.437   Cond. No.                         38.0\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "formula = '期末テスト ~ 小テスト + 睡眠時間'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.373829Z",
     "start_time": "2018-08-20T18:22:38.371457Z"
    }
   },
   "outputs": [],
   "source": [
    "eps_hat = np.array(result.resid)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 正規性の検定"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.378294Z",
     "start_time": "2018-08-20T18:22:38.374923Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.660"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.skew(eps_hat)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.382228Z",
     "start_time": "2018-08-20T18:22:38.379225Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.437"
      ]
     },
     "execution_count": 45,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.kurtosis(eps_hat, fisher=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ダービン・ワトソン比"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.386546Z",
     "start_time": "2018-08-20T18:22:38.383242Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.508"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.sum(np.diff(eps_hat, 1) ** 2) / np.sum(eps_hat ** 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-06-08T06:26:23.414638Z",
     "start_time": "2018-06-08T06:26:23.409582Z"
    }
   },
   "source": [
    "### 多重共線性"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.395450Z",
     "start_time": "2018-08-20T18:22:38.387540Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>小テスト</th>\n",
       "      <th>期末テスト</th>\n",
       "      <th>睡眠時間</th>\n",
       "      <th>通学方法</th>\n",
       "      <th>中テスト</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4.2</td>\n",
       "      <td>67</td>\n",
       "      <td>7.2</td>\n",
       "      <td>バス</td>\n",
       "      <td>8.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>7.2</td>\n",
       "      <td>71</td>\n",
       "      <td>7.9</td>\n",
       "      <td>自転車</td>\n",
       "      <td>14.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.0</td>\n",
       "      <td>19</td>\n",
       "      <td>5.3</td>\n",
       "      <td>バス</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3.0</td>\n",
       "      <td>35</td>\n",
       "      <td>6.8</td>\n",
       "      <td>徒歩</td>\n",
       "      <td>6.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>1.5</td>\n",
       "      <td>35</td>\n",
       "      <td>7.5</td>\n",
       "      <td>徒歩</td>\n",
       "      <td>3.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   小テスト  期末テスト  睡眠時間 通学方法  中テスト\n",
       "0   4.2     67   7.2   バス   8.4\n",
       "1   7.2     71   7.9  自転車  14.4\n",
       "2   0.0     19   5.3   バス   0.0\n",
       "3   3.0     35   6.8   徒歩   6.0\n",
       "4   1.5     35   7.5   徒歩   3.0"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df['中テスト'] = df['小テスト'] * 2\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-08-20T18:22:38.411290Z",
     "start_time": "2018-08-20T18:22:38.396466Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>          <td>期末テスト</td>      <th>  R-squared:         </th> <td>   0.676</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.658</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   37.61</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Tue, 21 Aug 2018</td> <th>  Prob (F-statistic):</th> <td>8.59e-06</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>03:22:38</td>     <th>  Log-Likelihood:    </th> <td> -76.325</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    20</td>      <th>  AIC:               </th> <td>   156.7</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    18</td>      <th>  BIC:               </th> <td>   158.6</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     1</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   23.6995</td> <td>    4.714</td> <td>    5.028</td> <td> 0.000</td> <td>   13.796</td> <td>   33.603</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>小テスト</th>      <td>    1.3107</td> <td>    0.214</td> <td>    6.133</td> <td> 0.000</td> <td>    0.862</td> <td>    1.760</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>中テスト</th>      <td>    2.6215</td> <td>    0.427</td> <td>    6.133</td> <td> 0.000</td> <td>    1.723</td> <td>    3.519</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 2.139</td> <th>  Durbin-Watson:     </th> <td>   1.478</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.343</td> <th>  Jarque-Bera (JB):  </th> <td>   1.773</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td> 0.670</td> <th>  Prob(JB):          </th> <td>   0.412</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.422</td> <th>  Cond. No.          </th> <td>1.09e+17</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                  期末テスト   R-squared:                       0.676\n",
       "Model:                            OLS   Adj. R-squared:                  0.658\n",
       "Method:                 Least Squares   F-statistic:                     37.61\n",
       "Date:                Tue, 21 Aug 2018   Prob (F-statistic):           8.59e-06\n",
       "Time:                        03:22:38   Log-Likelihood:                -76.325\n",
       "No. Observations:                  20   AIC:                             156.7\n",
       "Df Residuals:                      18   BIC:                             158.6\n",
       "Df Model:                           1                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     23.6995      4.714      5.028      0.000      13.796      33.603\n",
       "小テスト           1.3107      0.214      6.133      0.000       0.862       1.760\n",
       "中テスト           2.6215      0.427      6.133      0.000       1.723       3.519\n",
       "==============================================================================\n",
       "Omnibus:                        2.139   Durbin-Watson:                   1.478\n",
       "Prob(Omnibus):                  0.343   Jarque-Bera (JB):                1.773\n",
       "Skew:                           0.670   Prob(JB):                        0.412\n",
       "Kurtosis:                       2.422   Cond. No.                     1.09e+17\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "[2] The smallest eigenvalue is 1.65e-31. This might indicate that there are\n",
       "strong multicollinearity problems or that the design matrix is singular.\n",
       "\"\"\""
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "formula = '期末テスト ~ 小テスト + 中テスト'\n",
    "result = smf.ols(formula, df).fit()\n",
    "result.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  },
  "toc": {
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": "block",
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
